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ABSTRACT 



The problem of finding the shortest path between two 
points in the plane containing obstacles is considered. The 
set of such paths is uncountably infinite, making an 
exhaustive search impossible. This difficulty is overcome by 
reducing the size of the search space. The search is first 
restricted to a countably infinite set by focusing attention 
on the set of homotopy classes. By applying simple optimality 
principles, we obtain a finite list of such classes whose 
union contains the shortest path. This process of 
simplification is discussed in the thesis of CAPT Kevin D. 
Jenkins, U.S. Marine Corps. In this thesis we first discuss 
a computational investigation of two methods by which homotopy 
classes can be named. Next, a computational heuristic is 
presented that finds the lower bound for a path in a class. 
Finally, the true shortest path is found by searching these 
classes in order of increasing lower bound. One application 
of this study is in the area of robotic path planning. 
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THESIS DISCLAIMER 

The reader is cautioned that computer programs developed 
in this research may not have been exercised for all cases of 
interest. While every effort has been made, within the time 
available, to ensure that the programs are free of 
computational and logic errors, they cannot be considered 
validated. Any application of these programs without 



additional verification is at the risk of the user. 
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INTRODUCTION 



I . 

The industrial community relies today, in increasing 
measure, upon robots to handle a multitude of tasks. A need 
remains to develop robots, and fixed robot arms, with the 
ability to roam freely among obstacles. In this regard, the 
calculation of shortest paths is of obvious importance. 

A. THE PROBLEM 

As robot arms move, they must vary their configurations to 
position their end-effector or "hand". Once obstacles are 
introduced, freedom of movement of the manipulator may become 
drastically reduced. The problem considered here is one of 
finding the shortest path from a starting point to a 
destination point in a planar region containing obstacles. 
This paper addresses several parts of the solution to this 
probl em. 

B. OVERVIEW OF THE SOLUTION 

Between any two points in the plane, a path joining them 
may be chosen from an uncountably infinite set of 
alternatives. It is desirable to choose the shortest path 
from this set which, due to the presence of obstacles, may not 
be a straight line. An exhaustive search of the collection of 
possible paths is impossible, so another method must be 
developed . 
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The set of possible paths is partitioned into a collection 
of equivalence cl asses- -mutual 1 y exclusive subsets which 
collectively exhaust the partitioned set. This partition is 
induced when an equivalence relation is defined on the set of 
possible paths. This set of equivalence classes is countably 
infinite . 

Naming conventions are established which associate with 
each class a character string which allows the classes to be 
represented in a computer. 

Next, a finite list of candidate classes is produced to 
search. 

A heuristic is then applied to this finite list of classes 
to place them in an order which saves computational effort. 

The final step begins with this ordered list of candidate 
classes and solves the shortest path problem class by class in 
listed order. Savings of computational effort are realized 
when the class-by-class solution process can terminate without 
exhausting the ordered list of candidate classes. 

The key ideas in the above overview are now considered in 
somewhat greater detail. 

C . THE APPROACH 

1. The Topology 

Let P be the uncountably infinite set of all 
continuous, obstacle-avoiding paths from a starting point to 
a destination point, denoted a and z respectively. If p^ and 
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Pj are two paths in P, which have the same starting point and 



destination point, we say that p^ is homotopic to Pj if p^ can 
be mapped to Pj under a continuous function (with both 
endpoints fixed in place) without encroaching on any obstacles 
[Ref. l:p. 223]. Clearly, the homotopy relation is reflexive, 
symmetric and transitive and therefore defines an eguival ence 
relation [Ref. l:p. 223]. This relation induces a partition 
of P into a countably infinite collection of equivalence 
classes, known as homotopy classes [Ref. l:p. 223]. 

2. Establishing the Naming Convention 

In order to name homotopy classes, a reference frame 
is established to represent the topological relationships in 
the plane with obstacles. For a given path, p, a string of 
characters, R(p), is recorded which encodes information 
concerning the relationship of the path to the obstacles. 

Two algorithms are then presented which accept R(p) as 
input. These algorithms have the following important property 
which implicitly defines names for the homotopy classes: If 
P^ and Pj are coterminal paths and R(p^) and R(Pj) are input to 
either of the two algorithms, then the outputs are identical 
if and only if p^ is homotopic to pj [Ref. 2]. 

Chapter III presents a computational investigation 
which closely examines the methods employed by these two 
algorithms to produce the names for these homotopy classes. 
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3. Generation of Names for Candidate Equivalence Classes 
An edge-labeled graph is constructed which models the 

topological relationships within the region, where the nodes 
represent subregions induced by the reference frame. Nodes 
are connected if their corresponding subregions are adjacent. 
The names of the desired homotopy classes may be obtained by 
traversing this graph. This traversal produces a list of 
candidate classes which contain all paths of minimal length. 

The location of the origin affects the manner in which 
the plane is divided into wedges. In this regard, the graph 
which is created, given the location of the origin in one 
instance, may not be the same as the graph obtained if the 
origin is moved to a new point. 

4. Heuristic Ordering of Candidate Homotopy Classes 
For each class on the list, a lower bound on the 

length of its shortest representative path is constructed. 
The list of classes is then arranged into increasing order of 
these bounds. To obtain the bounds, a point is first chosen 
within each obstacle. A contraction deformation is then be 
applied to each obstacle to "shrink" the obstacle to that 
chosen point. 

Now some class is fixed and its shortest path is 
examined. As all obstacles are simultaneously contracted to 
their representing points, this shortest path has a limiting 
position which is polygonal . The length of this polygonal 
path is the lower bound which is associated with that class. 
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This length can be calculated from the class name and the 
representing points without explicitly defining any 
contraction deformation. This thesis presents this method to 
determine a bound for the shortest path in a given class. 

5. Class by Class Solution of the Shortest Path Problem 
In the final step of the solution, the classes on the 
ordered list of candidates are considered. The first class is 
removed from the list--that class with the smallest associated 
bound. The true shortest path in the named class is found by 
reversing the contraction previously applied to the obstacles, 
thereby transforming the polygonal path whose length provided 
the lower bound into the true shortest path. If the length of 
this path is smaller than the bound associated with the class 
on the top of the remaining list, the search is stopped. 
Otherwise, the first class is removed from that remaining list 
and the procedure is repeated. This procedure, which is 
contained in the thesis, continues until the condition 
specified above is met, and the shortest path has been found. 

D. Summary 

The solution to the problem of searching for the shortest 
path between two points in the plane with obstacles begins 
with consideration of a set of paths which is uncountably 
infinite. Through the homotopy relation, this set is reduced 
to a countably infinite set. The new set is further reduced 
to a finite list by modeling the region with a graph and 
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applying an optimality principle. This final list contains 
only those homotopy classes containing paths which are not 
sel f -crossing . From this list a shortest path is found. 

Computational effort is further reduced through the use of 
a heuristic which orders the list of candidate classes by 
increasing order of their lower bounds. The heuristic used 
also facilitates a methodical search of the classes while 
solving the shortest path problem. The use of this heuristic 
does not, however, imply that the solution is approximate. 
The solution to the problem will be exact using this method. 

E. THE CONTRIBUTION 

Chapter II of this thesis presents the methods used to 
establish a reference frame given the plane with obstacles and 
an algorithm which is used to generate the homotopy class name 
of a given path. 

Chapter III presents a computational investigation which 
addresses the question of whether or not two algorithms which 
are used to name homotopy classes actually provide the same 
results . 

Chapter IV introduces a method to find a lower bound for 
the shortest path in a class. 

Chapter V shows how to find the true shortest path in a 
class . 
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II. NAMING EQUIVALENCE CLASSES 



A. INTRODUCTION 

In this chapter, a labelling scheme is established to 
represent the equivalence classes. This notation is used 
throughout the paper to organize the computational search for 
the shortest path. 

An algorithm is presented that for a character 
representation of a path, p, names the path's respective 
equivalence class. It can be shown that, after being 

processed by the algorithm, different character 
representations of paths p^ and Pj yield the same output 
exactly when and only when the two paths are in the same 
equivalence class [Ref. 2]. Thus, the name of a class will be 
the string obtained by applying the algorithm to any 
representative in the class. 

The computation of class names depends on a reference 
frame which in turn depends on a collection of obstacles. 
Although more than one reference frame can be drawn for a 
particular collection of obstacles, the choice of a particular 
frame fixes the representation of all homotopy classes. 

Once a reference frame is developed, the name of the 
homotopy class for a path can be determined by a two-step 
procedure. First, a character string R(p) is calculated which 
encodes certain information about the path taken through the 
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obstacles. Second, this character string is accepted as input 
to an algorithm which then produces a name for the equivalence 
class to which the particular path belongs. This algorithm 
produces the class names in terms of the same alphabet used to 
create the initial character string. The results from the 
algorithm described. Algorithm 1, are used throughout this 
analysis as the class name associated with a given path. 

B. CONSTRUCTION OF REFERENCE FRAME 

Let b|j be an arbitrary point in obstacle B]^, k = 1,2, . . . ,n, 
where n is the number of obstacles in the region. A point c 
is chosen and a line drawn through each bj^, infinite in extent 
in each direction and having the following properties: there 
is an open disk, 6, centered at c such that 6nBjj=0 for all 
k=l,...,n, and the n lines connecting c with the points bj^ are 
distinct. Such a point c can always be found as the above two 
conditions are satisfied by any point in the planar region 
that is neither on an obstacle nor on the n(n-l)/2 lines 
determined by pairwise choices of distinct bj^. 

To draw a reference frame, the n lines are first 
constructed joining c with each h^. The line from c to bj^ is 
labeled as L|^ . Each line is then partitioned into two 
directed, semi-infinite rays and one finite length line 
segment. The ray directed from c in the direction away from 
b|j is called . The ray emanating from bj^ and away from c is 
denoted fl|j. The remaining line segment, [c,b|j], is also 
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denoted a|^ . The reference frame is the collection of line 
segments and rays so constructed, as illustrated in Figure 

2 . 1 . 




Figure 2.1 A Reference Frame With n=2 Obstacles 
C. RAW CHARACTER STRINGS 

A reference frame is fixed in a region T. Let a and z be 
points in T and let p be a directed path in T from a to z. 
Then the raw string of p, denoted R(p) , is defined to be the 
ordered sequence of characters obtained by following p from a 
to z and recording, in order, the names of the rays that are 
crossed . 

Two special cases must be addressed so as to make the 
above definition complete. First, in the case that p crosses 
no reference rays, we let R(p) = e, where e denotes the empty 
string. Second, if p crosses through c (simultaneously 
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crossing all O; ) the names of all will be recorded in order 

j J 

of increasing subscript. 

Thus, raw strings are of the form 

R(p) = x^xj. . . Xj, 

where each Xj belongs to the alphabet 

The character Xj (j = 1, 2,..., m) is the name of the 

reference ray crossed by p. Figure 2.2 shows a pair of paths 
connecting the two points a and z in the region with two 
obstacl es . 
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a 



1 



Figure 2.2 Paths a 2 B| and a ^02 



D. ALGORITHM 1 

Algorithm 1 accepts as input a raw string R(p) for any 
path p in T and any reference frame. The output is a string 
denoted C(R(p)) of characters also chosen from the alphabet A. 
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This output string C(R(p)) is the name of the homotopy class 
to which p belongs. 

The algorithm is presented in terms of two functions. The 
first is the sorting function o. Let S = Xj . . . Xjj be any string 
over A. If S contains a two character substring X]jX],^j = 
with i<j, then o(S) is the string which results by reversing 
the order of the leftmost such substring. Figure 2.3 depicts 
two such strings and makes clear that such paths are 
homotopic. If S contains no such two character substring, 
then o(S) = S. 



Z« 




Figure 2.3 Homotopic Alpha Strings a 2 aj and a ^02 

Repeated application of o sorts all substrings that 
consist entirely of into non-decreasing order of 

subscript. 2(S) is defined to be the string which results 
when this sorting is complete. 

The second function is the cancellation function %. If 
string S contains a two character substring with Xj^ = X|,^j 
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then x(S) is the string which results by removing the leftmost 
such two character substring. Otherwise x(S) = S. Such 
cancellation is intuitive, for a reference ray that is crossed 
twice consecutively is equivalent to one that is not crossed 
at all. Thus repeated application of x cancels all pairs of 
adjacent like occurrences of identical characters. Let X(S) 
be defined as the string which results when all such possible 
cancellation is complete. 

With these definitions complete. Algorithm 1 is given in 
Figure 2.4. 



begin 
k *- 0 
S,j - R(p) 

WHILE S|; not equal to X(2(S,j)) 
k - k+1 

Sj, - X(2(S,j.i)) 

END WHILE 
C(R(p)) - S|5 
END 



Figure 2.4 Algorithm 1 



The output from Algorithm 1, C(R(p)), is called the 
canonical representation of p and is the unique name for the 
homotopy class of p [Ref. 2]. That name is used throughout 
this paper. 
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III. 



A COMPUTATIONAL INVESTIGATION 



A. INTRODUCTION 

This chapter describes a computational investigation which 
was conducted as an aid to verifying that the homotopy class 
names produced by Algorithm 1 are equivalent to those produced 
by a known method of naming homotopy classes [Ref. 2]. 

A second algorithm. Algorithm 2, is introduced here to 
forge a link between the class names of Chapter II and the 
well-known fundamental group. Two processes are presented to 
exploit this link: If Algorithm 1 truly names the classes, 
then both procedures should produce the same output for every 
path tested. One million test cases were examined and no 
counter-examples were found. While this computational 
evidence is not a proof, it did provide initial support before 
the proof of the claim that Algorithm 1 names homotopy classes 
was found. [Ref. 2] 

B. ALGORITHM 2 

1. Fundamental Group 

The class name obtained from Algorithm 1 was written 
in terms of the alphabet A defined in Chapter II. We present 
an algorithm here where the name obtained is expressed in 
terms of the fundamental group of a topological space T. 
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The basic idea for producing a fundamental group for 
this problem is to regard the paths in T as elements of the 
group and path concatenation (*) as the group operation. 
There are two problems with considering the paths as the group 
elements. First, it is not necessarily possible to 
concatenate any two paths in T. In order to concatenate two 
paths using the operation (*), the first path must end at the 
point where the second begins. In this analysis the test 
paths all begin and end at the point a. This point is called 
a base point and concatenation is possible for every pair of 
paths . 

A second issue is that an inverse is not wel 1 defined 
for paths. If the identity is defined as the empty string, 
i.e., stay at the base point, then when a path and its inverse 
are concatenated we have traced out a path which is not equal 
to the identity. To avoid this problem, the homotopy classes 
of paths are considered. Then all paths that have an empty 
character string representation will be in the identity class. 
Also, any class followed by its inverse will equal the 
identity. In Figure 3.1, g^ represents the class aiCjf and 
(gj)’^ represents the class So, gj * (gj)’^ represents the 
class which reduces to the identity class by the rules 
of Algorithm 1. 

The group elements, therefore, are homotopy classes of 
loops around obstacles based at a common point and the group 
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operation is defined in terms of concatenation (*) of these 
classes . 




: Class 

g ^ ; class 

g ^ ^ g"' : P,a,a,p,= e 



Figure 3.1 A Homotopy Class and its Inverse 

Given a base point a, we let p be any loop beginning 
and ending at a. We let [p^] denote the homotopy class of P' . 
So the set of group elements is {[p^] such that p^ is a member 
of T, where p- is a loop based at a}. This set will be called 
G with the elements denoted gj. Figure 3.2 illustrates 
several paths and their homotopy classes given n=2 obstacles. 
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Then the group operation (*) can be defined by 

[P] * [q] = [pq]. 

It is important to note that the fundamental group is 
finitely generated. As Figure 3.3 illustrates, the generators 
are not unique. However, by fixing a set of generators, the 
class names become fixed. 
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a 



2 



Figure 3.3 Alternate Generators For Fundamental Group of the 

Space With Two Obstacles 

With the above information describing the fundamental 
group representation and a few necessary definitions. 
Algorithm 2 is given below. 
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2. Functions in Algorithm 

a. Side Array 

Consider the reference lines Lj^ constructed in 
Chapter II to be oriented rays with direction from c to b]j. It 
then becomes reasonable to discuss a 'right' and a 'left' side 
of those lines. A moving point is allowed to trace a path p 
from start point to destination point. Side(k) is a function 
which defines the side of L|^ on which the moving point lies. 
The output is either 'left' or 'right'. The output is never 
'on' because this routine is used only after a complete 
crossing of L]j takes place. A crossing is considered to be 
complete when the moving point leaves L|j to one side after 
having met Lj^ from the opposite side. 

b. Switch Function 

The switch function is defined by: 

{ 'right' if side (k) = 'left' 
switch (side(k)) = < 

I 'left' if side (k) = 'right' 

As the moving point traces p, each time a reference line Lj^ is 

crossed, the switch function will be applied to indicate on 

which side of the moving point lies. 

c. Index Function 

Let X]j be any character from the alphabet A 
representing the k^^ element in R(p). Let index(x|j) = j if Xj^ 
equals Oj or Bj. Thus, the k^*“ crossing completed by the moving 
point is a crossing of Lj . So if the moving point is on the 
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left side of Lj initially and this line is now crossed on 
either side of bj, then side(3) would equal 'right', 
d. The Algorithm 

Using the functions described above. Algorithm 2 
is given in Figure 3.4. The algorithm works in two phases. 
The first phase initializes the array side(k) for k = 1, . . . ,n. 
The second phase reads the raw string from left to right and 
adds an element to the fundamental group representation for 
each B crossing. When the end of the input string is reached, 
the output is a shorter string with one character 
corresponding to every B element in the original raw string. 
Each character represents a generator or its inverse. 



begin 

input R(p) = XjXj • . . 
for j = 1 , . . . , n 

If a is to right of reference line Lj then 
side (j) - right 

else 

side (j) •- left 
end if 
end for 
G - [a] 

for k = 1 , . . . , m 
i = index(xjj) 

side(i) *- switch(side(i) ) 
if X]^ = B^ for some r, then 

if side ( index ( X|j) ) = left, then 
G G * g^ 

else 

G - G * (g:)'! 
end if 
end if 
end for 
F(R(p)) - G 
end 

Figure 3.4 Algorithm 2 
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e. Fundamental Group Cancellation Function 

The cancellation rules in the fundamental group 
differ somewhat from those for the raw string. Although sets 
of generators are not unique, we can obtain a unique 
representation of each class with respect to a particular set 
of generators. For any given set of generators, every 

homotopy class can be represented as a product of these 
generators and their inverses. Even so, this representation 
is only unique after cancellation is applied. The 

cancellation rule follows. 

Let G = {gj, . . . , g^} be a set of generators of the 
fundamental group of T (base point a) and let Y = yjy 2 --*y[n be 
a representation of some homotopy class in terms of the gj and 
their inverses, i.e. for each i = 1 ,..., m, y^ = gj or (9][)'^ 
for some j,k = 1 ,..., n. 

The cancellation function K is defined as follows. 
If Y contains a two character substring YiYiti with y^ = 
then k(Y) is the string which results by removing both y^ and 
y’^j, in their leftmost occurrence, from Y. Finitely many 
repeated applications of K produce a string in which no 
further cancellation is possible. We define this string as 
K(Y). 
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C. THE COMPUTATIONAL INVESTIGATION 
1. The Approach 

A graphical representation of the two procedures is 
provided in Figure 3.5. The investigation accepts a random 
path p as input and generates two names for the homotopy class 
representing p as output. 



ALG 2 



F(C(R(p))) 



I 



R(p) 




k(F(R(p))) 

K(F(R(p))) 



Figure 3.5 Flow Chart of Two Algorithms 

The idea of the test is to determine for a given path 
its raw string character representation. For computational 
purposes, the paths considered are polygonal paths. The raw 
string is then used as input to both algorithms. 
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The raw string R(p) is first input to Algorithm 2. 
The output F(R(p)) from Algorithm 2 is then input to the 
cancellation routine and the fully cancelled string K(F(R(p))) 
is output. In parallel, R(p) is input to Algorithm 1 which 
produces the canonical string C(R(p)). Algorithm 2 is then 
applied to C(R(p)) to output F(C(R(p))). The results of these 
tests are then compared. In Figure 3.5 a flow chart 
illustrates this process. The program continues to compare 
the two outputs until a predetermined number of paths are 
checked . 

It can be shown that C(R(p)) is a unique class name if 
and only if F(C(R(p))) equals K(F(R(p))) for all p in 
T [Ref. 2]. Before the proof of this conjecture was obtained, 
the above procedures were programmed and successfully tested 
for one million different paths. 

2 . The Test 

The general algorithm used to test the two procedures 
represented in Figure 3.5 is easily followed and is shown in 
Figure 3.6. However, the attached program, which was used to 
test the model, is substantially more involved. This more 
involved program was written in an attempt to generate the 
most efficient programming code possible. 
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n = number of reference frames to be considered 
m = number of paths to be considered 
nobs = number of obstacles on each board 



begin 

for board = 1 , . . . , n 

create reference frame or board 
for path = 1 , . . . , n 

1. Create a polygonal path p with nseg 
segments 

2 . Form R( p) 

3. Form F(R(p)) 

4. Form K(F(R(p) ) ) 

5 . Form C(R(p) ) 

6. Form F(C(R(p) ) ) 

if K(F(R(p))) does not equal F(C{R{p))) then 
print (raw string, board, and path info 
for the counter example) 
end i f 
end for 
end for 

print (final seed) 
end 

Figure 3.6 Test Algorithm 



Two time saving techniques are employed in order to 
test the one million cases in 25 minutes. The first method 
involves the subroutines which perform cancellation. It is 
simple to code a program that scans the character 
representation repeatedly to find all possible cancellation in 
a string such as , which reduces to Bp However, in 
the enclosed code (Appendix A), pointers are inserted in the 
strings to mark the position where the first cancellation 
occurs. Then the newly adjacent characters are checked for 
possible cancellation. This method reduces the complexity of 
the test and is employed in both Algorithm 1 and Algorithm 2. 
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The second time-saving device was developed around the 
sorting of character strings. Instead of using a bubble sort 
method (order n complexity), a merge sort algorithm (order 
n(log(n)) complexity) is used. Table 1 shows the savings is 
realized immediately as opposed to there existing some cut off 
where n(ln(n)) becomes less expensive. 



TABLE 1 

COMPLEXITY COMPARISON 


n 




n(ln(n)) 


2 


4 


1 


5 


25 


8 


10 


100 


23 


100 


10,000 


460 


1,000 


1,000,000 


6,907 


1,000,000 


10E16 


13,815,510 



As can be seen in the attached code, this choice to 
save computer time requires a significant increase in program 
complexity and therefore in programmer time. The trade-off of 
computer time versus programmer time needs. to be considered in 
any similar analysis. 

The code used to implement the test is presented in 
Appendix A. 
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IV. DETERMINING A BOUND FOR A HOMOTOPY CLASS 



A. INTRODUCTION 

For each homotopy class on a list, a lower bound for the 
length of its shortest representative path is found. The list 
is obtained through a graph traversal process [Ref. 3]. The 
list of classes is then arranged in increasing order based 
upon these bounds. 

To obtain the bounds, we first fix a class and consider 
its shortest path, p. For example, consider the path in 
Figure 4.1 to be a string attached at z. The shortest path is 
the one obtained by pulling the string tight so that it lies 
against the obstacles and passes through "a". 





Figure 4.1 Class Representatives 



The topological deformation of "contraction" may then be 
applied to each obstacle to "shrink" the obstacle to an point 
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bjj which is chosen arbitrarily within the obstacle Bj^ [Ref. 4]. 
As seen in Figure 4.2, when all obstacles can be 
simultaneously contracted to bj^ . The shortest path p* found 
in this way has a limiting position which is polygonal. 




Figure 4.2 Class Representative to a Bound for 
the Shortest Representative 



The length of this polygonal path provides a lower bound that 
is associated with that class. The fact that this is a lower 
bound follows from the following argument. If we let D be the 
set of points through which paths may travel with full-sized 
obstacles and E be the set of available points in the space 
after the obstacles are "shrunk”, then D is a subset of E. 
Let n(X, E) be the set of all paths in equivalence class X 
which remains entirely in E. Since D is a subset of E, II(X,D) 
is a subset of II(X,E) for all X. We now define a function f 
that maps II(X,E) to a real number in R, with the rule that 
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f(p) = length of p. Since D is a subset of E, this equation 
foil ows : 



min f(p) min f{p) 
peO.{X,E) ^ peJl{x,D) 



{Figure 4.3) . 



This lower bound, minpjjj| gj f ( p ) , can be calculated from the 
class name and obstacle-representing points without explicitly 
defining any contraction deformation. 



f 




Figure 4.3 Mapping Related to Class Lower Bound 



Knowing that this polygonal path, p*, exists, we now give 
a method to find the minimum length path and compute its 
length. Effectively, we must connect the points bj^ in the 
correct order and compute the length of each line segment. 

The idea employed to find p* is: We are given a starting 
point "a" and a character string which represents the 
mandatory order of ray crossings. Each character of the 
string implies that we must cross a certain ray. We then 
compute the maximum number of rays (i.e. read the maximum 
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number of characters in the string) that can be crossed with 
a straight line. If z cannot be reached with a straight line 
we determine where it is best to make a turn. The turning 
point then becomes the new starting point and the process 
repeats until z is reached. To find the straight line 
segments of p* , the notion of a cone of directions in which 
the segments of a polygonal p* must lie is introduced. 

B, FINDING THE BOUND 

1. Establishing a Cone of Directions 

The cone of directions as depicted in Figure 4.4, is 
an open region--not including the boundaries-- formed by two 
rays based at a common point and forming an angle of less than 
71 radians. In this region we can draw any ray based at the 
point common to the boundary rays. This is the cone of 
directions . 




Figure 4.4 Cone of Directions 

2. Class Name Defines the Cone 

To create a cone we need a starting point "a" and a 
reference ray (Figure 4.5). We define an entire class name as 
Ojj for a path to cross the ray ctj^. First we draw a ray from 
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the starting point to the point representing the obstacle 



point b|j. In Figure 4.5, to cross ajj we must be to the right 
of the ray from a to b|j. The ray from "a" to bj^ is one 
boundary for the cone. Next we draw a ray parallel to the ray 
that we must cross, a^. Again in Figure 4.5, in order to 
cross with a straight line we must leave "a" at an angle 
that will allow us to cross the aj^ ray. Therefore, the ray 
parallel to forms the other boundary. The two boundaries 
form the cone of directions as any ray strictly lying between 
these two limits will cross ai^. 



b 

Starting 
point 

Figure 4.5 A Cone of Directions Toward 

3. Intersecting Multiple Cones 

Given a cone of directions for a single ray, we can 
then consider longer class name strings. For this we generate 
an intersection of multiple cones of direction. For example, 
in Figure 4.6, "a" is the starting point and the class is 
a2fl^i32 leading to the ending point z. 
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“l 

Figure 4.6 The Class a2fijfl2 

The cone of directions in which the first segment of 
p* must lie is determined by considering the first element of 
the class name. Therefore, the first cone must be in a 
direction which crosses the 02 reference ray. 

We draw a line from the starting point a to b2 and a 
second line from "a" parallel to 02 and in the same direction 
as 02 (Figure 4 . 7 ). This forms the cone of directions in 
which the first polygonal line segment of the shortest 
representative must lie. 




Figure 4.7 Cone of Directions Toward 0^2 
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Then we consider the second character and create 
the cone of directions that will insure crossing of B^ We 
intersect the two cones of direction. The cone of 
intersections is pictured in Figure 4.8. 




new cone 




Figure 4.8 



Two Intersecting Cones of Direction 



Continuing in this 
third character in 
cone . 



fashion, in Figure 4.9, we consider the 
the class, B 2 and create the appropriate 
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empty 

intersection 



Figure 4.9 An Empty Intersection of Directions 



We see that the intersection of these cones is empty. 
To enter this new cone we are forced out of the cone in which 
we know the first segment must lie. Therefore, we choose a 
path which crosses all rays that were named ( 0:2 and 13^) and 
turns on some bj^. We determine the correct bj^ below. 

4. An Empty Intersection 

As can be seen in Figure 4.9, the new cone of 
directions lies to the right of the previous intersected cone 
of directions. This implies that the path we require will 
bend to the right. The peg--point b]^ on which the path bends- 
-is found on the right hand boundary of the old cone. In 
fact, as shown in Figure 4.9, the obstacle centered at b^ is 
the peg that the path bends on in the class a2D^B2- 
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We now know that one segment of the bounding path is 
the line from a to . In the next step we use bj as the new 
starting point. Since 02 have been crossed, we use B2 
to create the cone of directions in Figure 4 . 10 . 




z does not 
lie in the cone 



Now, since there are no more characters to read, we 
check to see if z lies within the most recently established 
cone of directions. In Figure 4.10 we see that z does not lie 
in the cone so the bounding path must bend again. 

Since z lies to the right of the cone, the bounding 
path bends on b2 which is on the right hand boundary of the 
cone of directions. Hence, the final two segments of the 
bounding path are the line segments between bj and b2 and the 
line segment between b2 and z. 

So given the class <12^1^2 obtain the bounding path 
pictured in Figure 4.11 on the contracted obstacles. 
Computing the length of this path is now straight forward. 
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The next section generalizes these steps and shows how 



to find the path to z. 




Figure 4.11 The Bounding Path for the Class a2Bjfl2 



C. THE ALGORITHM 



To find a lower bound for a path in a given class, the 
following algorithm is applied: 



t( j) : = 



peg(i ) 
cone ( t ( 



U : = 




L,R : = 



Lg , Rj, : - 



an array of points with t(l)=a, and t(m)=z; 
t(i+l) through t(m-l) are pegs visited along the 
bounding polygonal path. 

:= the obstacle-representing point associated 
with Xj 

j),x^):= the cone of directions from a starting 
point t(j) in which we can cross the ray 
x^ with a straight line from t 
represents the cone associated with the most 
recently read character 

records the most recent cone of intersections 
represents the cone of directions created by one 
character; V is assigned to V^j^j and is replaced by 
the intersection of U and 

represent pointers into the character string X 
and is associated with the obstacle on the left 
of right boundary of the cone of directions. If 
no boundary lies on a ray, then L or R is 
assigned zero 

associated with the obstacle on the left of right 
boundary of the cone U that is intersected with 
Vgjj. If one ray of U lies inside V^ij, then either 
L or R is replaced by Lg or Rg . It both rays of 
U lie inside then L and R are replaced by Lg 

and Rg respectively. 
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BEGIN 

t(l) <— a; i <-- 1 
procedure START 
read character x{ 

initialize current cone V <— cone(t, x[) 

IF peg( i ) on left boundary of V THEN 
L <- i; R <- 0 
ELSE 

R <“ i; L <“ 0 
END IF 

procedure CONTINUE 
i <- i+1; read character 
U <— cone(t, xj) 

IF peg(xi) on left boundary of V THEN 
Lu <- i; Ru 0 
ELSE 

Ru <-- i; Lu <-- 0 
END IF 

^old ^old ^Id ^ 



V <- U n Void 

IF left boundary of U falls within Void THEN; L <— Lu; END IF 
IF right boundary of U falls within Void THEN; R <-- Ru; END IF 
procedure BRANCH 
IF V is empty THEN 

IF U is to left of Void THEN 
j<-j + l; t (j)<--peg(Loid) 

ELSE ( U is to right of Void ) 
j<--j + l; t (]■)<-- peg(Rold) 

END IF 

go to procedure START 
ELSE (V is not empty.) 

IF i < m (not all characters have been read) THEN Go to procedure CONTINUE. 
ELSE (There are no more characters to read.) 



IF z lies in current cone V THEN j <— j + 1 
ELSE 

IF U is to left of Void THEN 


; tO')<--z 


; STOP 


j <-- j + 1; t(j) <- peg(Loid); j <-- j + 1 
ELSE ( U is to right of Void ) 


; t(j) <" z; STOP 


j <- j + 1; t(j) <- peg(Roid); j <-- j + 1 
END IF 


; tO)<--z 


; STOP 



END IF 
END IF 
END IF 
END 
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V. FINDING THE SHORTEST PATH 



A. INTRODUCTION 

In order to find the true shortest path from a to z we 
consider the classes of paths on the ordered list of candidate 
paths. We proceed by removing the first class from the list. 
The true shortest path in the named class is found. If the 
length of this path is smaller than the bound associated with 
the class on the top of the remaining ordered list, the search 
is stopped. Otherwise, the first class is removed from the 
remaining list and the above procedure is repeated. This 
continues until the true shortest path is found. 

To find the shortest path within a given class, we reverse 
the contraction applied to the obstacles performed earlier 
when we found the lower bound for the classes (see Figure 
5.1). For each polygonal path, whose length provided a lower 
bound, we transform the segment into the true shortest path in 
that class. Accordingly, the process of calculating the 
shortest path in a class begins with the polygonal path used 
to lower bound it and a description of the obstacle 
geometries . 
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cx. 




1 



Figure 5.1 Reversing the Obstacle Contraction 



To this point, no assumptions have been made about the 
shapes of the obstacles. In this chapter, with little loss of 
generality and for computational purposes which will be 
discussed later, all obstacles are assumed to be circles 
centered on bj, . We could, of course, assume other shapes, 
like polygons or ellipses. 

B. FINDING LINE SEGMENTS 

As can be seen in Figure 5.2 the line segment from a to Bj 
is tangent to B^. Also, the line segment from Bj to the 
obstacle B 2 is tangent to both obstacles. In general there 
are two tangent lines from a point to a circle and there are 
four distinct lines tangent between two circles (see Figure 
5.2). In this section we see how to find these tangent lines 
and we show how to pick the correct tangent line to represent 
the true shortest path in the class. 
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Figure 5.2 Possible Tangent Lines 

Assuming that the shortest path from a to z in a given 
class of paths is not a straight line, we use another 
application of the cone of directions idea to determine on 
which side of an obstacle to turn. That is equivalent to 
picking the correct tangent line. We use the last cone of 
intersections formed just before the intersection of cones 
becomes empty to make this determination. 

In Figure 5.3 we see that the path is required to bend 
around obstacle Bj . If a tangent line is drawn from a to the 
left side of obstacle 1 and we then travel on an arc around 
the obstacle until we can proceed directly to z, we would have 
the shortest path in the class and could easily find its 
length. If not then we iterate the process. 
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Figure 5.3 Path Bend in Class 



The following method is used to find the tangent lines to 
a circular obstacle when given the starting point, the 
coordinates of the center of the obstacle, and the radius of 
the obstacle. 

In Figure 5.4, the starting point, a, is the point (u,v). 
The center coordinates and the radius of the peg are the point 
(h,k) and r, respectively. 

y 



Figure 




5.4 Point to Obstacle Tangent 



Lines 
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with this information we determine the equations of the 



circle 


representing the 


obstacle and an 


arbitrary 


line. 


We 


solve 


these equations 


simultaneous 1 y 


for the 


points 


of 



intersection (Figure 5.5). The specific value that we want is 
the point of tangency. 





Figure 5.5 The Intersections of a Circle and a Line 



So, given h, k, and r and assuming, without loss of 
generality, that b=u=v=0, we find the following two equations: 
(x-h)’ + (y-k)^ = r* 

y = mx (where m is the slope of the line) 

By expanding the first equation and substituting y=mx we 
obtain : 

(l+m")x" + (2mk-2h)x + (h"+k"-r") = 0 . 

By the quadratic formula we obtain: 

^ ^ -2mk+2h± [ {2mk-2h) ^-4 ( 1 +m^) (h^+k^-r^) ] 

2 (1+CJ2) 



In Equation (1), if the discriminant is zero, we have found 
the value where the points of intersection coincide. So 
setting the discriminant equal to zero yields: 
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(-2A-2icm)2-4 (l+/n2) (/3*+Jc«-r2) =0 
(A+icm)2-(m8+i) (ii*+ic*-r2) =0 
(A*/ri*+2AJan+tf) -(/ri*+l) (/i^+ic^-r*) =0 
(k!^-l^-kf^+r^)nF+{2hk) m+(r^-l^) =0 



from which m 
A second 



can be determined. 

application of the quadratic formula 



-2hk± [4/22/c 2-4 (r2-A2) (r^-A*) ] 
2 (r*-i3*) 






-hk±r [72*+ic*-r*] ^ 
r^-t^ 



i=l,2 . 



produces : 



In general, with an arbitrary u and v 



in,- = 



_ - (h-u) {k-v) ±r [ (h-u) (k-b) ^-r^] 



Vi 



r2-(A-u)2 



i=l,2 



The only time this technique does not appear to work is 
when the tangent line found is a vertical line (m is undefined 
because there is a zero in the denominator). This problem can 
be easily handled because if m is undefined we know the 
equation for a vertical line to be x equal to some constant. 

The value for m is then substituted into Equation (1) to 
find the , i = 1,2, coordinate and the value of x^ is 
substituted into y^ = m^x^ to find y^. By knowing on which side 
of the obstacle the path must turn, we can pick the tangent 
line with the appropriate slope. A method of picking the 
correct line will be discussed in the next section. 
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We follow the arc along the circle to the point of 
tangency of the line between the current obstacle and z or to 
the point of tangency to the next obstacle. Tangent lines 
between two obstacles can be determined by a similar method as 
the one just shown. A method of finding the four tangent 
lines between two circular obstacles follows. 

C. DETERMINING THE CORRECT TANGENT LINES 
1. Point to Obstacle 

We have shown how to find the two tangent lines from 
a point to a circle. We present a method to determine which 
tangent line to pick. 

The obstacle on which we turn lies on a boundary of 
the cone of directions. Thus the cone itself lies either to 
the left or to the right of the ray through the obstacle 
center. If the cone lies on the left, use the tangent of the 
left is taken. If the cone lies to the right, then the 
tangent on the right is taken. The case with the cone on the 
left is illustrated in Figure 5.6. 




Figure 5.6 Use the Left Tangent for Path Bending to the Right 
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2. Obstacle to Obstacle 



When considering two circular obstacles, there are 
four distinct lines that can be drawn tangent to both 
obstacles. In this section we show how to pick the tangent 
line that corresponds to the respective shortest path in the 
given class. 

To find the four distinct tangent lines between two 
obstacles we first create a line which is tangent to one of 
the circles. We then move this tangent line along the circle 
boundary until we discover the four points at which the line 
is tangent to both obstacles (Figure 5.7). 




Figure 5.7 Rotating 4> to Find Tangent Lines Between 
Two Circles 
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Without loss of generality, in Figure 5.7 we fix one 
obstacle center at the origin. We pick an arbitrary point 
(x',y’), on the obstacle represented as (r*cos(<|>), 
r*sin(<|))). We construct a line through (x',y') tangent to Bp 
The slope m of this line is: 

y-r*sin (<|>) 
x-r*cos{<j>) 



This slope is also equal to -cot(4>). Hence, 



m 



= -cot (4.) 

x-r*cos4> 



( 2 ) 



The equation of the circle with center (h,k) is 
ix-h)^+(y+k)^ = i?2 

Rearranging Equation (2) we obtain the following: 

y - r*sin(4>) = (-cot (4>) ) (x - r* (cos (4>) ) 
y= (-cot (4>) )x + r*cos (4>) *cot (4>) + r*sin(<|>) 



Equation (3) is in slope-intercept form and the y intercept b 
is : 
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b = r*cos (4>) *cot (4)) + r*sin(<|>) 



“ (i^ * ) 

_ / cos2<l) + sin*4> \ 

■ [ sin<|> J 

b = r*csc(4>) 

Thus both m and b are functions of <|). The general line that 
we use to intersect both circles is 

y - in(^)x + jb(4>) 

This equation for the line and for the obstacle centered at 
(h,k) are solved simultaneously. After steps similar to those 
used in the point to obstacle example we obtain the following 
f ormul a : 

(/7^+1)a® + [-2h+2m(b-k)]x + [b^+ {b-k)*-B^] = 0 . 

Again by applying the quadratic formula and setting the 
discriminant equal to zero we find 

A[-h-m(b-k)]^ - [b^+(b-k)^-m =0 . 

Substituting the function values for m(<j>) and b(<|)) yields 

[-A- (cot<|>) (r*cso(j>-ic) ]*-(cot*4>+l) [li^+ (r*csc^-k)^-F^] =0 (4) 
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This equation has four zeros for values of (J) in [0,27t) 



as can be seen in Figure 5.8. 







Figure 5.8 Four Zeros of Equation (4) 

Since <() is the only parameter in Equation (4) we can 
use the four values found to label the four tangent lines. 
The line associated with the lowest value of 4> is labeled A, 
the second line associated with the second lowest (J) value is 
labeled B and the last two lines are labeled with C and D, 
respective! y . 
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Assuming the polygonal path bends on each peg on the 
path, there are essentially two possible arrangements when 
considering a polygonal path that negotiates two obstacles. 
These possibilities are shown in Figure 5.10. The paths can 
be flipped about the horizontal axis to produce two equivalent 
paths . 





Figure 5.10 Two Possible Polygonal Paths 



46 



In order to discover which situation exists, we extend 
the middle line segment and use it as a reference in Figure 
5.10. Since the path bends at each peg, the two end segments 
must lie on one side or the other of the reference line. If 
they both lie on the same side, then we have situation a 
above. If the two segments lie on opposite sides, then we 
have situation b. 

Given the two possible arrangements in Figure 5.10, we 
know the corresponding shortest paths in the class look like 
those in Figure 5.11. 




Figure 5.11 True Shortest Paths 
We have shown how to pick the correct tangent line 
from a point to a circle. This is applied to find the first 
and last segments of this polygonal path. 

To find the correct tangent between two obstacles, 
following the situation on the left in Figure 5.11 we must 
pick tangent line B, described earlier. If we have the other 
situation in Figure 5.11, we must pick tangent line A. 

If the polygonal path does not bend on each peg, then 
there are two possible configurations as seen in Figure 5.12. 
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Figure 5.12 Polygonal Paths With Collinear Segments 



To this point, we only used the cone of directions 
idea to find the correct tangent from a point to an obstacle. 
We use the cone of directions information to pick the correct 
tangent lines too. 

As seen in Figure 5.13 we consider the class B^a 2 . 




The tangent line from a to is found in the same 
manner as discussed earlier. However, since the path of the 
bound does not bend, we cannot choose the correct tangent line 
between and B 2 in the same way as before. 

Instead, we construct a reference line that connects 
the centers of the two obstacles. Then we consider on which 
side of this line the cones of directions lie (Figure 5.14). 



If they lie on different sides of the reference line then we 
pick tangent line A or D. Otherwise, we pick tangent lines B 
or C . 




Figure 5.14 Cone of Direction to Help Pick Tangent Lines 

D. TWO POTENTIAL PROBLEMS TO CHECK 

Given a path in a class like the one in Figure 5,15, we 
see the associated bounding path when we shrink the obstacles 
to points. However, there are two separate problems that may 
arise when the obstacles are restored to their original shape: 
1) The prospective path may pass through one or more 
obstacles. 2) When the obstacles are restored to their full 
size the path created can switch sides of one or more 
obstacles. These possibilities can be seen in Figure 5,16. 
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Figure 5.15 A Class and Its Bounding Path 
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the bounding 
path 




a 






2 




Figure 5.16 A Bounding Path and Its Restored Path 

The path determined by the cone of directions technique 
passes through obstacles and Bj, illustrating the first 
problem. By switching sides of Bj, the path is no longer in 
the same equivalence class, illustrating the second problem. 
Actually, the path should look like the one in Figure 5.17 
which avoids all obstacles and remains in the original class. 
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Figure 5.17 Corrected Shortest Path 

The correct path around the obstacles can be determined by 
using a distance function. Then, methods discussed in Chapter 
IV are employed to negotiate the new impeding obstacles. In 
Figure 5.18, segment Pj*--which is one segment of the bounding 
path--joins pegs bj^ and bj^^, and lies along the line through 
these pegs. We let a circular obstacle Bj have center (Uj,Vj) 
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and radius . Pegs bj^ and have cartesian coordinates 

(x^,y-) and respectively. We let d(*,*) denote 

Euclidean distance. 

The distance from p^* to Bj--which is the perpendicular 
distance from p* to the nearest point in Bj--is given by 

d(p*,B-) if the center of B.: lies in the band denoted in Figure 

5 J J 

5.18 whose borders are perpendicular to p^* and which passes 
through bj^ and bj^^j , and by min{d(b]j , Bj) , d( bj^^j , B') } if not. We 
can then check whether the center of Bj is in the band if we 
define vectors 

= ([Xi^-x^], [y^ti-yi]) 

and 

U'j = ([u-X'], [V'-y^]) 

which can be seen in Figure 5.19 [Ref. 5]. 




Figure 5.18 Band Created By One Segment of a Polygonal Path 
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Then the center of B. is in the band if and only if X; • ja; , ^>0 
and X' • y^j<0, where (•) denotes the inner product. We write 
the distance function corresponding to Pj* passing through Bj 
as 

f d{p*B,hj) ; if Xj • j>0 and X^ • 

min d{by.,Bj) ,d{bj^^^,Bj) ; otherwise 



These distances are easily computed with 



and 



dipl.Bj) 









[Ref. 5]. 




Figure 5.19 



Vectors 



to Help Determine 



d(Pj*,Bp 
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With this new information, appropriate tangent lines must 
be found in order to avoid newly-violated obstacles and to 
keep the path in the correct class. These new tangent lines 
may also encroach obstacles. So the process of checking all 
obstacles is continued until no obstacles are violated. 

To find the length of the shortest path, we sum the 
straight line segments and all of the arc lengths 
corresponding to the distances traveled around each circular 
obstacle . 
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VI. CONCLUSIONS AND RECOMMENDATIONS 



A. CONCLUSIONS 

Determination of the shortest path between two points in 
a plane containing obstacles has applications in many fields, 
particularly that of robotic path planning. This thesis 
addresses several key issues belonging to a larger plan for 
the solution of the shortest path problem. This plan may be 
summarized as follows: 

1. Partition the set of paths into homotopy classes. 

2. Associate algebraic "names" with each of these classes. 

3. Produce a finite list of feasible classes. 

4. For each class on the above list, calculate an 
inexpensive lower bound on the length of its shortest 
path . 

5. Arrange the list in increasing order of these lower 
bounds . 

6. Considering the classes in listed order, calculate the 
shortest path in the class. Calculation terminates when 
some class is encountered whose shortest path has length 
shorter than the lower bound associated with the 
succeeding class. 

The approach outlined above offers two primary 

computational advantages. First, note that the collection of 
possible paths from a given point of origin to a given 
destination is uncountably infinite. By organizing the search 
around homotopy classes and by taking advantage of certain 
topological relationships in the region [Ref. 3], the search 
of an uncountably large set is replaced by a search of a 
finite list. 
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Second, the search of this finite list is ordered 
according to lower bounds which estimate the lengths of the 
associated shortest paths. This ordering allows the search to 
terminate without exhausting the list. The bounds associated 
with each class are inexpensive to calculate--relative to the 
cost of finding the class' shortest path. This results in 
further savings of computational effort. 

The contributions of this thesis concentrate in areas 2, 
4, and 6 of the above list. 

As regards the association of names with homotopy classes 
(area 2) we present a computational investigation which was 
conducted ancillary to proof activity in this area. In 
summary, we associate a name with each homotopy class as 
follows: by imposing a certain labeled reference frame on the 
region in question; by encoding information about the 
relationship of a given path with the frame; and by employing 
algorithms presented here. In this way, we reduce the encoded 
information to a string of characters which is unique to each 
class. The investigation presented here uses a presentation 
of the well-known fundamental group as a standard of 
comparison to verify that the character strings used in our 
approach are in fact class names. (This proposition has 
subsequently been proved [Ref. 2]). The significance of the 
structure of the names which we employ is that they 
incorporate information which is useful both for the 
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calculation of class bounds and for calculation of class 



shortest paths. 

For the calculation of class bounds (area 4) we represent 
each obstacle by a point chosen within it. This point is used 
to construct the reference frame which is imposed on the 
region. Restricting our attention to such summary information 
allows class bounds to be calculated cheaply. 

The final step (area 6) begins with an ordered list of 
classes. In this paper we present a method by which shortest 
paths of a class may be calculated when the obstacles are 
circular. This is the only step in which any generality is 
lost, and is potentially fruitful subject for further study. 

B. RECOMMENDATIONS FOR FURTHER STUDY 

The field of robotics is one that is gaining much interest 
in mathematics and other fields of study. Further analysis of 
the reference frame and its properties may aid future 
researchers in answering other questions concerning the 
shortest path problem. 

One aspect of the problem that can be studied further 
relates to areas 4 and 6 above. Programming the two steps of 
finding the lower bound for a class and subsequently 
determining the shortest path in a class would be quite 
useful . 

Another computational problem relates to area 2. It is 
believed that further study and refinement of the code could 
result in a program which is even more efficient. 
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Another possible subject of future research is that of 
considering non-circular obstacles (area 6). It should be 
noted that any set of obstacles can be covered with circles. 
With this approximation the procedures presented in this 
thesis may be applied. 

Finally, this thesis introduces only part of the solution 
to the shortest path problem. The thesis by CAPT Kevin D. 
Jenkins, U.S. Marine Corps [Ref. 3], presents the remaining 
portions of the solution (area 3), and should be studied in 
conjunction with this paper. 
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APPENDIX A. FORTRAN PROGRAM: !>£ OOMHJTATIONAL INVESTIGATION 



PROGRAM CLASNAME FORTRAN 

Q *4c:»:4:*********4c**4c***4:***4:***4:*****4c*****4:******4c*:4c#*:t::t:**:|c:t;#**4::tc3^:^:C##:t::tc 

C THIS PROGRAM RANDOMLY GENERATES OBSTACLES AND A POLYGONAL PATH 
C THROUGH THOSE OBSTACLES TO DETERMINE WHETHER OR NOT THE PROCEDURES 
C USED BY COMPETING ALGORITHMS PRODUCE THE SAME FUNDAMENTAL GROUP 
C REPRESENTATION FOR THAT PATH. 

Q *:^4:t:^4:4:**4:**********4:4:*****4c*****#**:tc**:tc4:****4:***4:**4uMc**************4: 

C THIS PORTION OF THE PROGRAM SERVES AS THE MAIN DRIVER WHICH RECEIVES 
C THE PARAMETERS DEFINING THE REGION AND INITIALIZES THE LINK LIST 
C ARRAYS WHICH WILL BE USED TO REPRESENT THE POLYGONAL PATH. 

INPUT: INITIAL SEED FOR THE RANDOM NUMBER GENERATOR, NUMBER OF 

OBSTACLE CONFIGURATIONS AND PATHS TO BE TESTED, NUMBER OF 
OBSTACLES IN THE PLANE AND NUMBER OF SEGMENTS IN EACH 
POLYGONAL PATH 

OUTPUT: FINAL RANDOM NUMBER GENERATOR SEED AND MESSAGES INDICATING 

ANY PATHS WHICH PRODUCE DIFFERENT FUNDAMENTAL GROUP 
REPRESENTATIONS 

REAL*8 BX(IOOO), BY(IOOO), X(1000), Y(1000), DSEED 
INTEGER NOBS, NSEGS, NUMPTS, N, M, HEAD(IOOO) 

+ NEXTOOOO), PRED(IOOO), BOARD, PATH 

EXTERNAL QGUBFS 

DSEED = 123457.0 
M = 1000 
N = 1000 
NOBS = 20 
NSEGS = 5 
NUMPTS = NSEGS + 1 
PRINT*, ■ INPUT SEED’, DSEED 

DO 1 BOARD = 1, N 

CALL BOARDS (DSEED, BX, BY, NOBS) 

DO 2 PATH = 1, M 

CALL IN IT (HEAD, NEXT, PRED, PSEED, DSEED ) 

CALL PATHS (DSEED, X, Y, NSEGS) 

CALL TEST ( HEAD , NEXT , PRED , BX , BY , NOBS , X , Y , NUMPTS , BSEED , PSEED ) 
2 CONTINUE 
1 CONTINUE 

PRINT*,' ' 

PRINT*, 'FINAL SEED', DSEED 

STOP 
END 
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SUBROUTINE INIT ( HEAD, NEXT, PRED , PSEED, DSEED) 

Q t:^Pt**t**t***************t***t*******itf*****t*tt****^***it:**t***t******** 

C THIS SUBROUTINE INITIALIZES THE DOUBLE LINK LIST ARRAYS THAT WILL 
C REPRESENT THE PATH. 

C INPUT: DSEED 

C OUTPUT: HEAD, NEXT, AND PRED ARRAYS SET TO ZERO, AND SEED FOR THE 
C R.N.G. PRIOR TO CONSTRUCTION OF THE PATH 

INTEGER HEAD (1000), NEXT (1000), PRED (1000) 

REAL*8 PSEED, DSEED 

PSEED = DSEED 

DO 3 K = 1,1000 
HEAD(K) = 0 
NEXT(K) = 0 
PRED(K) = 0 
3 CONTINUE 

RETURN 

END 



SUBROUTINE BOARDS (DSEED, BX, BY, NOBS) 

C 

C THIS SUBROUTINE USES A PSEUDO RANDOM NUMBER GENERATOR TO CREATE THE 
C COORDINATES OF EACH OBSTACLE ON THE BOARD, SCALING ALL COORDINATES 
C TO BE IN THE INTERVAL (-1,1). 

C INPUT: NUMBER OF OBSTACLES IN THE REGION AND A SEED FOR THE R.N.G. 

C OUTPUT: OBSTACLE COORDINATES 

REAL*8 BX(IOOO), BY (1000), DSEED 
EXTERNAL GGUBFS 

DO 1 1=1, NOBS 

BX(I) = 2. * GQUBFS(DSEED) - 1. 

1 BY(I) = 2. * GGUBFS(DSEED) - 1. 

RETURN 

END 
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SUBROUTINE PATHS (DSEED, X, Y, NSEGS) 

C RANDOMLY GENERATES THE X AND Y COORDINATES OF THE VERTICES FOR THE 
C POLYGONAL PATH, SCALING ALL COORDINATES TO BE IN THE INTERVAL (-1,1). 
C THIS SUBROUTINE ALSO ENSURES THAT THE PATH IS A CLOSED LOOP BY 
C ASSIGNING THE START/FINISH POINTS THE SAME COORDINATES. 

C INPUT: NUMBER OF PATH SEGMENTS AND A SEED FOR THE R.N.G. 

C OUTPUT: COORDINATES OF VERTICES ALONG THE POLYGONAL PATH 

REAL*8 X(1000), Y(1000), DSEED 
EXTERNAL GGUBFS 

DO 1 1=1, NSEGS 

X(l) = 2. ♦ QGUBFS(DSEED) - 1. 

1 YO) = 2. * QGUBFS(DSEED) - 1. 

X(NSEGS+1) = X(1) 

Y(NSEGS+1) = Y(1) 

RETURN 

END 



SUBROUTINE TEST (HEAD, NEXT, PRED, BX, BY, NOBS, X, Y, NUMPTS, 

+ BSEED, PSEED) 

C GENERATES THE RAW STRING OF CHARACTERS REPRESENTING A PATH AND 
C DETERMINES ITS FUNDAMENTAL GROUP REPRESENTATION USING ALGORITHMS 1 
C AND 2. RESULTS OF THESE COMPETING ALGORITHMS ARE THEN COMPARED FOR 
C DIFFERENCES. 

INPUT: A POLYGONAL PATH AND COORDINATES OF POINTS REPRESENTING 

OBSTACLES 

C OUTPUT: COMPARED RESULTS OF ALGORITHMS 

INTEGER HEAD(IOOO), NEXT(IOOO), PRED(IOOO), FR(IOOO), NR(IOOO), 

+ PR(IOOO), FL(IOOO), NL(IOOO), PL(IOOO) 

REALMS BX(IOOO), BY(IOOO), X(1000), Y(1000) 

DATA FR/ 1000*0/ 

DATA NR/ 1000*0/ 

DATA PR/ 1000*0/ 

DATA FL/ 1000*0/ 

DATA NL/ 1000*0/ 

DATA PL/ 1000*0/ 
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CALL RAWSTR(BX, BY, NOBS, X, Y, NUMPTS, HEAD, NEXT, PRED, NELEMS) 
CALL ALG2 (NOBS, NELEMS, HEAD, NEXT, X, Y, BX, BY, FR, NR, PR) 
CALL CANALG2 (FR, NR, PR) 

CALL ALG1 (HEAD, NEXT, PRED, NELEMS) 

CALL ALG2 (NOBS, NELEMS, HEAD, NEXT, X, Y, BX, BY, FL, NL, PL) 
CALL CTREX (FL, NL, PL, FR, NR, PR, BSEED, PSEED) 

RETURN 

END 



SUBROUTINE RAWSTR(BX, BY, NOBS, X, Y, NUMPTS, HEAD, NEXT, PRED, 

+ NELEMS) 

Q 

C PRODUCES THE RAW STRING OF CHARACTERS WHICH REPRESENTS THE PATH, 

C PAYING CLOSE ATTENTION TO THE ORDERING OF THE CHARACTERS WHERE 
C NECESSARY. THE STRING IS CONSTRUCTED BY IDENTIFYING THOSE OBSTACLE 
C REFERENCE LINES WHICH ARE CROSSED AS EACH PATH SEGMENT IS TRAVERSED 
C I N ORDER . 

INPUT: ORDERED LIST OF VERTICES REPRESENTING THE PATH AND A POINT 

REPRESENTING EACH OBSTACLE IN THE REGION 

OUTPUT: A RAW STRING OF ALPHAS AND BETAS CONTAINED IN AN ARRAY 

NAMED HEAD() AND ITS PARALLEL ARRAYS NEXT() AND PRED() WHICH 
PRODUCE THE DOUBLE LINKED LIST 

IMPLICIT REAL*8 (A-H,0-Z) 

INTEGER HEAD, NEXT, SEG, FSTSTR, HDINDX, SEGEND, FIRST, START, 

+ LENGTH, PRED 

LOGICAL ALLALF 

DIMENSION X(1000), Y(1000), BX(IOOO), BY(IOOO), A1(1000), 

+ B1(1000), A2(1000), B2(1000), D2(1000), DIST(IOOO), 

+ FIRST(IOOO), HEAD(IOOO), NEXT(IOOO), PRED(IOOO) 

DATA D I ST/ 1000*0.0/ 

DATA FIRST/ 1000*0/ 

CALL SETUP(NUMPTS,N0BS,CX,CY,X,Y,BX,BY,A1 ,B1 , A2, B2,D2) 

HDINDX = 2 
NSEGS = NUMPTS - 1 
LOLD = 1 
HEAD(LOLD) = 0 
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DO 7 SEG = 1, NSEGS 
SEGEND = SEG + 1 
FSTSTR = HDINDX 
LSTSTR = HDINDX 
ALLALF = .TRUE. 

LENGTH = 0 

DO 6 L I NE = 1 , NOBS 

CHECK1 = (A1(LINE)*X(SEG)+B1(LINE)*Y(SEG))* 

+ (A1(LINE)*X(SEGEND)+B1(LINE)*Y(SEGEND)) 

IF (CHECK1.lt. 0) THEN 
LENGTH = LENGTH + 1 
CHECK2 = (A2(SEG)*BX(LINE)+B2(SEG)* 

+ BY(LINE)+D2(SEG))*D2(SEG) 

IF (CHECK2.lt. 0) THEN 
HEAD(HDINDX) = -LINE 
LSTSTR = HDINDX 
HDINDX = HDINDX + 1 

ELSE 

CALL CASES1 (Al, B1, A2, B2, D2, SEG, LINE, XI NT, 
+ YINT) 

DISTC = XINT**2 + YINT**2 

DISTB = (XINT-BX(LINE))**2 + (YINT-BY(L INE))**2 
IF (DISTB.lt. DISTC) THEN 
HEAD(HDINDX) = LINE 
LSTSTR = HDINDX 
HDINDX = HDINDX + 1 
ALLALF = .FALSE. 

ELSE 

HEAD(HDINDX) = -LINE 
LSTSTR = HDINDX 
HDINDX = HDINDX + 1 
END IF 
END IF 
END IF 

6 CONTINUE 

IF (LENGTH. NE.O) THEN 
IF (ALLALF) THEN 

CALL ALPHAS (NEXT, LOLD, FSTSTR, LSTSTR, SEG, START, 

+ HEAD, NSEGS, PRED) 

ELSE 

CALL ORDER(SEG, LINE, HDINDX, HEAD, NEXT, Al, A2, 

+ B1, B2, D2, FSTSTR, LSTSTR, X, Y, START, 

+ LENGTH, FIRST, LOLD, NSEGS, PRED) 

END IF 
END IF 

7 CONTINUE 

CALL OOUNTR( START, NEXT, HEAD, NELEMS) 

RETURN 

END 
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SUBROUTINE SETUP ( NUMPTS , NOBS, CX, CY, X, Y, BX, BY, Al, B1, A2, 

+ B2, D2) 

Q 

C FOR EACH OBSTACLE, THIS ROUTINE DETERMINES THE COEFFICIENTS OF THE 
C EQUATION FOR THE REFERENCE LINE FROM THE OBSTACLE TO THE ORIGIN. IN 

C ADDITION, IT CALCULATES THE COEFFICIENTS OF THE LINE REPRESENTING 
C EACH SEGMENT OF THE POLYGONAL PATH. 

C INPUT: NUMPTS, NOBS, COORDINATES OF VERTICES ALONG POLYGONAL PATH 

C AND COORDINATES OF THOSE POINTS WHICH REPRESENT EACH 

C OBSTACLE 

C OUTPUT: COEFFICIENTS OF LINEAR EQUATIONS REPRESENTING PATH SEGMENTS 

C AND REFERENCE LINES FOR EACH OBSTACLE 

IMPLICIT REAL*8 (A-H, 0-Z) 

DIMENSION X(1000), Y(1000), BX(IOOO), BY(IOOO), A1(1000), 

+ B1(1000), A2(1000), B2(1000), D2(1000) 

DO 4 I = 1, NOBS 
A1(l) = BY(I) 

Bl(l) = -BX(I) 

4 CONTINUE 

CX = 0.0 

CY = 0.0 

I = 1 

DO 5 J = 2, NUMPTS 

A2(l) = Y(J) - Y(l) 

B2(l) = X(l) - X(J) 

D2(l) = (Y(I)*X(J))-(X(I)*Y(J)) 

I = 1 + 1 

5 CONTINUE 

RETURN 

END 
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SUBROUTINE CASES1(A1, B1, A2, B2, D2, SEG, LINE, XINT, YINT) 

0 



C THIS SUBROUTINE DETERMINES THE COORDINATES FOR THE POINT OF 
C INTERSECTION OF A GIVEN PATH SEGMENT AND A GIVEN OBSTACLE REFERENCE 
C LINE. NUMERICAL STABILITY OF CALCULATIONS REQUIRES THE MATHEMATICAL 
C OPERATIONS TO BE SEPARATED INTO CASES, THE BEST CASE BEING USED FOR 
C EACH PARTICULAR SITUATION. 



C INPUT: COEFFICIENTS OF LINEAR EQUATIONS FOR PATH SEGMENT AND 

C REFERENCE LINE TO BE EVALUATED 



C OUTPUT: POINT OF INTERSECTION OF THE PATH SEGMENT AND THE OBSTACLE 

C REFERENCE LINE 

IMPLICIT REAL*8 (A-H,0-Z) 

INTEGER SEG 

DIMENSION A1(1000), B1(1000), A2(1000), B2(1000), D2(1000), A(2,2) 



A(1,1) = Al(LINE) 

A(1,2) = BI(LINE) 

A(2,1) = A2(SEG) 

A(2,2) = B2(SEG) 

BIGEST = 0.0 

DO 1 L = 1,2 

DO 1 K = 1,2 

TEST = DABS(A(K,L)) 

IF (TEST. GT. BIGEST) THEN 
BIGEST = TEST 
KBIG = K 
LBIG = L 
END IF 

1 CONTINUE 



IF (KBIG.EQ.1) THEN 

IF (LBIG.EQ.1) THEN 

YINT = -D2(SEG)/(B2(SEG)-B1(LINE)*A2(SEG)/A1(LINE)) 
XINT = -B1(LINE)*YINT/A1(LINE) 

RETURN 

ELSE 

XINT = -D2(SEG)/(A2(SEG)-A1(LINE)*B2(SEG)/B1(LINE)) 
YINT = -A1(LINE)*XINT/B1(LINE) 

RETURN 
END IF 

ELSE 

IF (LBIG.EQ.1) THEN 

YINT = (D2(SEG)*A1(LINE)/A2(SEG))/ 

+ (B1(LINE)-B2(SEG)*A1(LINE)/A2(SEG)) 

XINT = (-D2(SEG)-B2(SEG)*YINT)/A2(SEG) 

RETURN 
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GO GO 



ELSE 



XINT = (D2(SEG)*B1(LINE)/B2(SEG))/ 

+ (A1(LINE)-A2(SEG)*B1(LINE)/B2(SEG)) 

YINT = (-D2(SEG)-XINT*A2(SEG))/B2(SEG) 
RETURN 
END IF 
END IF 

END 



SUBROUTINE ALPHAS (NEXT, LOLD, FSTSTR, LSTSTR, SEG, START, HEAD, 
+ NSEGS, PRED) 

C GIVEN THAT A SEGMENT OF THE PATH CROSSES ONLY ALPHA RAYS WHEN 
C TRAVERSED, THESE CROSSINGS ARE SIMPLY INSERTED INTO THE LINK LIST 
C ARRAYS IN THE ORDER IN WHICH THEY WERE DETECTED, I.E., SMALLEST TO 
C LARGEST IN ABSOLUTE VALUE. 

INPUT: STRING REPRESENTING THE REFERENCE LINES CROSSED BY THE 

CURRENT PATH SEGMENT, GIVEN THAT ALL CROSSINGS ARE ALPHAS 

OUTPUT: UPDATED NEXT() AND PRED() ARRAYS WHICH CONTAIN THE STRING 

REPRESENTING THE MOST CURRENT PATH SEGMENT 

INTE(XR NEXT(IOOO), FSTSTR, SEG, START, HEAD(IOOO), PRED(IOOO) 

NEXT(LOLD) = FSTSTR 
PRED (FSTSTR) = LOLD 
LAST = LSTSTR - 1 

IF ( LAST. GT. FSTSTR) THEN 
DO 12 I = FSTSTR, LAST 
NEXT(I) =1+1 
PREDO + 1) = I 
12 CONTINUE 

ELSE IF (LAST. EQ. FSTSTR) THEN 
NEXT( FSTSTR) = LSTSTR 
PRED (LSTSTR) = FSTSTR 
END IF 

IF(SEG.EQ.I) START = LOLD 

LOLD = LSTSTR 
RETURN 
END 
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SUBROUTINE ORDER(SEG, LINE, HDINDX, HEAD, NEXT, A1, A2, B1, B2, 

+ D2, FSTSTR, LSTSTR, X, Y, START, LENGTH, FIRST, 

+ LOLD, NSEGS, PRED) 

C GIVEN THAT A SEGMENT OF THE PATH CROSSES ONE OR MORE BETAS, THIS 
C SUBROUTINE DETERMINES THE ACTUAL ORDER OF CROSSING WHEN TRAVERSING 
C THE SEGMENT FROM BEGINNING TO END. THIS IS DONE BY FIRST DETERMINING 
C THE DISTANCE FROM THE SEGMENT START POINT TO THE POINT OF 
0 INTERSECTION OF EACH CROSSED OBSTACLE REFERENCE LINE. THESE 
0 DISTANCES ARE THEN SORTED FROM SMALLEST TO LARGEST AND CROSSINGS 
C ARE UPDATED IN THE LINK LIST ACCORDINGLY. 

C INPUT: HEADO AND ALL COORDINATES REQUIRED TO DETERMINE THE 

C DISTANCES FROM INITIAL VERTEX OF PATH SEGMENT TO EACH OF 

C THE CROSSED REFERENCE LINES 

C OUTPUT: HEADO, NEXT() AND PRED() ARRAYS CONTAINING THE RAW STRING 

C WHICH ACCURATELY REPRESENTS THE PATH ALONG THE CURRENT 

C SEGMENT 

INTEGER SEG, LINE, HEAD(IOOO), NEXT(IOOO), FSTSTR, HDINDX, 

+ START, FIRST(IOOO), PRED(IOOO), F 

REALMS XINTER(IOOO), Yl NTER( 1000) , AI(IOOO), B1(1000), 

+ A2(1000), B2(1000), D2(1000), DIST(IOOO), X(1000), Y(1000) 

DO 8 J = FSTSTR, LSTSTR 
LINE = ABS(HEAD(J)) 

CALL CASES2 (A1, B1, A2, B2, D2, SEG, LINE, XINTER, YINTER) 
DIST(J) = (XINTER(LINE)-X(SEG))**2+(YINTER(LINE)-Y(SEG))**2 
DIST(J) = -DIST(J) 

NEXT(J) = J+1 
PRED(J+1) = J 
8 CONTINUE 
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NEXT (LSTSTR) = 0 
PRED(LSTSTR+1) = LSTSTR 

CALL MERG2(DIST, NEXT, PRED, FSTSTR, F) 

NEXT(LOLD) = F 
PRED(F) = LOLD 
IF(SEG.EQ.I) START = LOLD 
LOLD = NEXT (LOLD) 

IF(NEXT(LOLD).NE.O) THEN 
LOLD = NEXT (LOLD) 

GO TO 17 
END IF 



RETURN 

END 
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GO GO 



SUBROUTINE CASES2(A1, B1, A2, B2, D2, SEG, LINE, XINTER, YINTER) 

Q 

C DETERMINES COORDINATES FOR THE POINT OF INTERSECTION OF THE PATH 
C SEGMENT AND EACH OF THE REFERENCE LINES IT CROSSES. AGAIN, IN ORDER 
C TO MAINTAIN NUMERICAL STABILITY, CALCULATIONS ARE MADE USING THE 
C MOST APPROPRIATE CLOSED FORM EQUATION. 

INPUT: COEFFICIENTS OF LINEAR EQUATIONS FOR CURRENT PATH SEGMENT 

AND ALL OBSTACLE REFERENCE LINES WHICH IT CROSSES 

OUTPUT: COORDINATES FOR POINTS OF INTERSECTION OF EACH OBSTACLE 

REFERENCE LINE WITH THE PATH SEGMENT 

IMPLICIT REAL*8 (A-H, 0-Z) 

INTEGER SEG 

DIMENSION A1(1000), B1(1000), A2(1000), B2(1000), D2(1000), 

+ XINTER(IOOO), YINTER(IOOO), A(2,2) 



A(1 


.1) 


= 


AI(LINE) 


A(1 


.2) 


z 


BI(LINE) 


A(2 


.1) 


z 


A2(SEG) 


A(2 


,2) 


= 


B2(SEG) 


BIGEST 


= 


0.0 


DO 


1 L 


z 


1,2 




DO 


1 


K = 1,2 



TEST = DABS(A(K,D) 

IF (TEST.GT.BIGEST) THEN 
BIGEST = TEST 
KBIG = K 
LBIG = L 
END IF 



1 CONTINUE 



IF (KBIG.EQ.1) THEN 

IF (LBIG.EQ.1) THEN 

Yl NTER(L I NE)= -D2(SEG)/(B2(SEG)-B1 (L INE)*A2(SEG)/A1 (L I NE) ) 
XINTER(LINE)= -B1 (L I NE)*Y I NTER(L I NE)/A1 (L I NE) 

RETURN 

XI NTER(L I NE)= -D2(SEG)/(A2(SEG)-A1 (L I NE)*B2(SEG)/B1 (LI NE)) 
Yl NTER(L I NE)= -Al (L I NE)*X I NTER(L I NE)/B1 (L I NE) 



ELSE 



ELSE 



RETURN 
END IF 

IF (LBIG.EQ.1) THEN 

YINTER(LINE) = (D2(SEG)*A1(L I NE)/A2(SEG) )/ 

(B1 (LI NE)-B2(SEG)*A1 (L I NE)/A2(SEG) ) 
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XINTER(LINE) = (-D2(SEG)-B2(SEG)*YINTER(LINE))/A2(SEG) 
RETURN 

ELSE 

XINTER(LINE) = (D2(SEG)*B1 (L INE)/B2(SEG) )/ 

+ (A1 (L I NE)-A2(SEG)*B1 (L I NE)/B2(SEG) ) 

YINTER(LINE) = (-D2(SEG)-XINTER(LINE)*A2(SEG))/B2(SEG) 
RETURN 
END IF 
END IF 
END 



SUBROUTINE MERG2 (DIST, NEXT, PRED, FSTSTR, F) 

C THIS SUBROUTINE SORTS A SUBSTRING OF ALL POSTIVE INTEGERS INTO 
C INCREASING ORDER. 

C INPUT: A DOUBLE LINK LIST CONSISTING OF DIST, NEXT, AND PRED ARRAYS 

OUTPUT: A DOUBLE LINK LIST WITH ALL ENTRIES PLACED IN THE ORDER 
IN WHICH THEIR RESPECTIVE REFERENCE LINES WERE CROSSED 

IMPLICIT INTEGER(A-Z) 

REAL*8 DIST (1000) 

LOGICAL DONE 

DIMENSION NEXT (1000), PRED(IOOO) 

DONE = .FALSE. 

P = FSTSTR 

1 F = P 

CALL SORT2 (F, PREDF, P, DIST, NEXT, DONE, PRED) 

IF(DONE) RETURN 
GOTO 1 
END 
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SUBROUTINE SORT2(F, PREDF, P, DIST, NEXT, DONE, PRED) 

3|c3|c sjc sjc sjc sjc sjc sjc s|C 3|Cs|C 3|C3|C 3|C3^9|( s|C sjCSjC sjCsjC 3|C3jCs|C 3|C s|C3|C 3|o|C sjcsjc sjcsjc sjC 3jC3|C 3|C 3|C3)C3jC s|C sjcsjc 3^3|(2^3|C 3^3|C sjcsjc s|C sjCSjC sjCSjC 3|C3|C s|C2^?^3|C3|C sjCSjC 3|C3|C3^3^3^ 

C MARCHES DOWN A SUBSTRING OF POSITIVE INTEGERS AND ONLY SORTS IF 
C ELEMENTS ARE NOT IN INCREASING ORDER. IF A NUMBER NEEDS TO BE 
C PLACED HIGHER IN THE LIST SUBROUTINE 'PUT' IS CALLED TO DO SO 

INPUT: POINTER F INTO DOUBLE LINK LIST ARRAYS DIST, NEXT, AND PRED 

TO INDICATE THE BEGINNING OF A SUBSTING OF POSITIVE INTEGERS, 
AND PREDF 

OUTPUT: P IS A POINTER, DIST(P) IS THE LAST POSITIVE INTEGER IN THE 
SUBSTRING THAT IS BEGUN BY DIST(F) 

IMPLICIT INTEGER (A-Z) 

REAL*8 DIST(IOOO) 

LOGICAL DONE 

DIMENSION NEXT(IOOO), PRED(IOOO) 

TAIL = F 

1 NTAIL = NEXT(TAIL) 

IF(NTAIL.EQ.O) THEN 
DONE = .TRUE. 

RETURN 
END IF 

IF(DIST(NTAIL).LE.DIST(TAIL)) THEN 
TAIL = NTAIL 
GOTO 1 
END IF 

CALL PUT2(F, TAIL, NTAIL, PREDF, DIST, NEXT, PRED) 

GOTO 1 
END 
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SUBROUTINE PUT2 (F, TAIL, NTAIL, PREDF, DIST, NEXT, PRED) 



C REARRANGES POINTERS TO PLACE POSITIVE INTEGERS IN INCREASING ORDER. 



C INPUT: F 

C TAIL 

C NTAIL 

C 
C 

C PREDF 



START OF POSITIVE SUBSTRING 
END OF CURRENTLY SORTED PORTION OF SUBSTRING 
NEXT(TAIL) - POINTER TO THE SUCCESSOR OF TAIL IN THE 
ITEM TO BE INCORPORATED INTO THE SORTED PORTION OF 
THE LIST 

THE PREDECESSOR OF F IN LINKED LIST 



OUTPUT: SOME POINTERS IN NEXT AND PRED ARE CHANGED TO PUT DIST(NTAIL) 
IN ITS PROPER PLACE IN THE SORTED PORTION ON THE LIST 



IMPLICIT INTEGER (A-Z) 

REALMS DIST (1000) 

DIMENSION NEXT(IOOO), PRED(IOOO) 

IF (DIST(NTAIL).GE.DIST(F)) THEN 
TEMP = NEXT (NTAIL) 

NEXT(NTAIL) = F 
PRED(F) = NTAIL 
NEXT(TAIL) = TEMP 
PRED(TEMP) = TAIL 
F = NTAIL 
RETURN 
END IF 

CALL WALK2(F, TAIL, NTAIL, DIST, NEXT, PRED) 

RETURN 

END 
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SUBROUTINE WALK2(F, TAIL, NTAIL, DIST, NEXT, PRED) 

Q Hp:^:;iiipt^t*******t*************************************************^***t* 

C THIS SUBROUTINE WALKS DOWN THE LINKED LIST AND PLACES THE DIST(NTAIL) 

C IN THE CORRECT POSITION IN THE DIST ARRAY. 

C INPUT: F, TAIL, NTAIL ARE USED AS ABOVE 

C OUTPUT: ALTERS POINTERS IN NEXT AND PRED ARRAYS TO PLACE DIST(NTAIL) 

C AFTER DIST(F) AND BEFORE DIST(TAIL) IN DOUBLE LINK LIST 

IMPLICIT INTEGER (A-Z) 

REAL*8 DIST(IOOO), HNTAIL 
DIMENSION NEXT(IOOO), PRED(IOOO) 

I = F 

NEXT I = NEXT(I) 

HNTAIL = DIST(NTAIL) 

1 IF(HNTAIL.QE.DIST(NEXTI )) THEN 
NEXT(I) = NTAIL 
PRED(NTAIL) = I 
NNTAIL = NEXT(NTAIL) 

NEXT(NTAIL) = NEXT I 
PRED(NEXTI) = NTAIL 
NEXT(TAIL) = NNTAIL 
PRED( NNTAIL) = TAIL 
RETURN 
END IF 

I = NEXTI 
NEXT I = NEXT(I) 

GOTO 1 

END 
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SUBROUTINE OOUNTR ( START , NEXT, HEAD, NELEMS) 

C THIS SUBROUTINE COUNTS THE NUMBER OF ELEMENTS I N A Ql VEN STR I NG OF 
C CHARACTERS 

C INPUT: HEADO AND NEXT() ARRAYS FOR THE STRING OF CHARACTERS 

C OUTPUT: NUMBER OF ELEMENTS IN THE STRING 

INTEGER PTR, NEXT (1000), HEAD (1000), START 

PTR = NEXT (START) 

NELEMS = 0 

20 IF (PTR.NE.O) THEN 

NELEMS = NELEMS + 1 
PTR = NEXT (PTR) 

GO TO 20 
END IF 

RETURN 

END 
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SUBROUTINE ALG2(NOBS, NELEMS, HEAD, NEXT, X, Y, BX, BY, FR, NR, 

+ PR) 

0 ****+**:):******* 4 :**:|:*************************************************** 

C GENERATES THE WELL-KNCW FUNDAMENTAL GROUP REPRESENTATION OF AN 
0 EQUIVALENCE CLASS. 

INPUT: RAW OR CANONICAL STRING, WITH # OBSTACLES, # ELEMENTS IN 

STRING, COORDINATES OF A AND OBSTACLES (BK) 

OUTPUT: 1) IF INPIJT IS RAW STRING, THEN POSSIBLY UNREDUCED 

FUNDAMENTAL GROUP WILL RESULT (IF CANCELLATION OF LIKE 
POSITIVE NUMBERS COULD HAVE OCCURRED IN THE RAW STRING). 

2) IF INPUT IS CANONICAL STRING, THEN THE RESULTING 
FUNDAMENTAL GROUP WILL BE IN REDUCED FORM. 

REAL*8 XA, YA, XB, YB, X(1000), Y(1000), BX(IOOO), BY(IOOO) 
INTEGER NOBS, NELEMS, HEAD(IOOO), NEXT(IOOO), PRED(IOOO), 

+ FR(IOOO), NR(IOOO), PR(IOOO), S(1000), FUNDGP( 1000) , 

+ START 

LOGICAL RIGHT(IOOO), RITE 

LENGTH = 1 
START = NEXT(1) 

XA = X(1) 

YA = Y(1) 

DO 1 K = 1, NOBS 
XB = BX(K) 

YB = BY(K) 

1 RIGHT(K) = RITE(XA, YA, XB, YB) 

NELEM = NELEMS + 1 
DO 10 M = 2, NELEM 

S(M) = HEAD(START) 

J = IABS(S(M)) 

RIGHT(J) = .NOT. RIGHT(J) 

IF(S(M).GT.O) THEN 

LENGTH = LENGTH + 1 
IF(RIGHT(J)) THEN 

FUNDGP( LENGTH) = J 

ELSE 

FUNDGP (LENGTH) = -J 

END IF 
END IF 

START = NEXT (START) 

IF(START.EQ.O) GOTO 11 
10 CONTINUE 
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11 FR(1) = 0 
NR(1) = 2 

DO 2 1 =2, LENGTH 
FR(I) = FUNDGP(I) 
NR(I) = I + 1 
PR(I) = I - 1 
2 CONTINUE 

NR (LENGTH) = 0 

RETURN 

END 



LOGICAL FUNCTION RITE(XA, YA, XB, YB) 

C THIS FUNCTION DETERMINES WHICH SIDE OF A GIVEN DIRECTED LINE ANY 
C POINT LIES. 

C INPUT: TWO POINTS THAT DETERMINE THE LINE 

OUTPUT: LOGICAL VARIABLE THAT IS TRUE IF A POINT LIES TO THE RIGHT 
AND FALSE IF A POINT LIES TO THE LEFT 

REAL*8 XB, YB, XA, YA, SIGNA 

RITE = .TRUE. 

SIGNA = -((YB*XA)-(XB*YA)) 

IF(SIGNA.LT.O)THEN 
RITE = .FALSE. 

RETURN 
END IF 

IF(SIGNA.EQ.O) THEN 

PRINT*, 'THE POINT A LIES ON THE LINE LK PROGRAM STOPS' 

STOP 
END IF 

RETURN 
END 



76 



o o 



SUBROUTINE CANALG2(HEAD, NEXT, PRED) 

C THIS SUBROUTINE TAKES AN UNREDUCED FUNDAMENTAL GROUP REPRESENTATION 
C OF A GIVEN CLASS AND CANCELS A GENERATOR IF IT IS ADJACENT TO ITS 
C I NVERSE . 

C INPUT: HEAD, NEXT, AND PRED ARRAYS 

OUTPUT; HEAD, NEXT, AND PRED ARRAYS WITH NEXT AND PRED REARRANGED 
TO SKIP AROUND CANCELLED ELEMENTS 

IMPLICIT INTEGER (A-Z) 

DIMENSION HEAD(IOOO), NEXT(IOOO), PRED(IOOO) 

START = 1 
PTR1 = START 
PTR2 = START 
PTR3 = NEXT(PTR2) 

10 IF (PTR3.NE.0) THEN 

IF (HEAD(PTR2).EQ.-(HEAD(PTR3))) THEN 
NEXT(PTRI) = NEXT(PTR3) 

PRED ( NEXT (PTR3)) = PTRl 
IF(NEXT(PTR1).EQ.O) RETURN 
PTR2 = PTRl 
PTRl = PRED(PTRI) 

PTR3 = NEXT(PTR2) 

ELSE 

PTRl = PTR2 
PTR2 = PTR3 
PTR3 = NEXT(PTR2) 

END IF 
GO TO 10 
END IF 
RETURN 
END 
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SUBROUTINE ALG1(HEAD, NEXT, PRED, NELEMS) 

C THIS SUBROUTINE TAKES A GIVEN RAW STRING OF CHARACTERS REPRESENTING 
C A PATH AND PRODUCES THE CANONICAL FORM OF THAT STRING. THIS IS DONE 
C BY FIRST ORDERING ALL OF THE ALPHA SUBSTRINGS FROM SMALLEST TO 
C LARGEST IN ABSOLUTE VALUE. NEXT, CANCELLATION IS PERFORMED TO 
C ELIMINATE ALL LIKE PAIRS OF ADJACENT ELEMENTS FROM THE STRING. 

INPUT: RAW STRING IN FORM OF DOUBLE LINKED LIST WITH HEAD(), NEXT() 

AND PREDO ARRAYS 

C OUTPUT: CANONICAL FORM OF THE RAW STRING 

IMPLICIT INTEGER (A-Z) 

DIMENSION HEAD(IOOO), NEXT(IOOO), PRED(IOOO) 

START = 1 

CALL MERG1(HEAD, NEXT, PRED) 

CALL CANCEL (START, HEAD, NEXT, PRED, NELEMS) 

RETURN 

END 
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SUBROUTINE MERG1 (HEAD, NEXT, PRED) 

Q 

C THIS SUBROUTINE DOES THE INITIAL SORT OF THE ALPHA SUBSTRINGS IN THE 
C RAWSTRING (INCREASING IN ABSULUTE VALUE). 

INPUT: HEAD, NEXT, PRED ARRAYS REPRESENTING A DOUBLE LINKED LIST 

OF THE RAWSTRING WITH THE ALPHA SUBSTRINGS UNORDERED 

OUTPUT: POINTERS STORED IN THE ARRAYS NEXT AND PRED ARE ALTERED SO 
THAT EACH SUBSTRING OF THE STORED LIST WHICH CONSISTS 
ENTIRELY OF NEGATIVE INTEGERS IS SORTED INTO NON- I NCREAS I NG 
ORDER, WHILE SUBSTRINGS OF POSITIVE INTEGERS ARE LEFT 
UNALTERED. 

IMPLICIT INTEGER(A-Z) 

LOGICAL DONE 

DIMENSION HEAD(IOOO), NEXT(IOOO), PRED(IOOO) 

DONE = .FALSE. 

P = 1 

1 CALL FRONT (P, F, PREDF, HEAD, NEXT, DONE) 

IF (DONE) RETURN 

CALL SORT (F, PREDF, P, HEAD, NEXT, DONE, PRED) 

IF (DONE) RETURN 
GOTO 1 
END 
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SUBROUTINE FRONT(P, F, PREDF, HEAD, NEXT, DONE) 

C ♦♦♦♦♦♦♦♦♦♦♦*****************************************5<c***!M:************ 

C FINDS THE BEGINNING OF NEGATIVE INTEGER STRINGS (ALPHA STRING) 

C INPUT: POINTER P INTO LINKED LIST 

OUTPUT: F IS A POINTER INTO TTHE HEAD ARRAY. POINTS TO FIRST NEGATIVE 
ENTRY WHICH OCCURS STRICTLY AFTER HEAD(P). PREDF IS POINTER 
SUCH THAT HEAD(PREDF) IS THE PREDESSOR OF HEAD(F). 

IMPLICIT INTEGER(A-Z) 

LOGICAL DONE 

DIMENSION HEAD(IOOO), NEXT(IOOO) 

F = P 

1 PREDF = F 
F = NEXT(F) 

IF (F.EQ.O) THEN 
DONE = .TRUE. 

RETURN 
END IF 

IF (HEAD(F).LT.O) RETURN 
GOTO 1 
END 
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SUBROUTINE SORT(F, PREDF, P, HEAD, NEXT, DONE, PRED) 

Q 

C MARCHES DOWN A SUBSTRING OF NEGATIVE INTEGERS AND ONLY SORTS IF 
C ELEMENTS ARE IN INCREASING ORDER. IF A NUMBER NEEDS TO BE 
C PLACED HIGHER IN THE LIST SUBROUTINE ’PUT’ IS CALLED TO DO SO 

INPUT: POINTER F INTO DOUBLE LINK LIST ARRAYS HEAD, NEXT, AND PRED 

TO INDICATE THE BEGINNING OF A SUBSTRING OF NEGATIVE 
INTEGERS; AND PREDF 

OUTPUT: P IS A POINTER, FEAD(P) IS THE LAST NEGATIVE INTEGER IN THE 
SUBSTRING THAT IS BEGUN BY HEAD(F) 

IMPLICIT INTEGER (A-Z) 

LOGICAL DONE 

DIMENSION HEAD(IOOO), NEXT(IOOO), PRED(IOOO) 

TAIL = F 

1 NTAIL = NEXT(TAIL) 

IF(NTAIL.EQ.O) THEN 
DONE = .TRUE. 

RETURN 
END IF 

IF(HEAD(NTAIL).GT.O) THEN 
P = TAIL 
RETURN 
END IF 

IF(HEAD(NTAIL).LE.HEAD(TAIL)) THEN 
TAIL = NTAIL 
GOTO 1 
END IF 

CALL PUT(F, TAIL, NTAIL, PREDF, HEAD, NEXT, PRED) 

GOTO 1 
END 
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SUBROUTINE PUT (F, TAIL, NTAIL, PREDF, HEAD, NEXT, PRED) 

Q ***********************************:4:*******:^c************************^* 

C SUBROUTINE THAT REARRANGES POINTERS TO PLACE A NEGATIVE INTEGER IN 
C NON- INCREASING ORDER. 



C INPUT: 

C 

C 

C 

C 

C 



F - START OF NEGATIVE SUBSTRING 
TAIL - END OF CURRENTLY SORTED PORTION OF SUBSTRING 
NTAIL = NEXT(TAIL) - POINTER TO THE SUCCESSOR OF TAIL IN THE 
ITEM TO BE INCORPORATED INTO THE SORTED PORTION OF 
THE LIST. 

PREDF - THE PREDECESSOR OF F IN LINKED LIST 



OUTPUT: SOME POINTERS IN NEXT AND PRED ARE CHANGED TO PUT HEAD( NTAIL) 
IN ITS PROPER PLACE IN THE SORTED PORTION ON THE LIST 



IMPLICIT INTEGER (A-Z) 

DIMENSION HEAD(IOOO), NEXT(IOOO), PRED(IOOO) 

IF (HEAD(NTAIL).GE.HEAD(F)) THEN 
NEXT(PREDF) = NTAIL 
PRED (NTAIL) = PREDF 
TEMP = NEXT(NTAIL) 

NEXT(NTAIL) = F 
PRED(F) = NTAIL 
NEXT(TAIL) = TEMP 
PRED (TEMP) = TAIL 
F = NTAIL 
RETURN 

END IF 

CALL WALK(F, TAIL, NTAIL, HEAD, NEXT, PRED) 

RETURN 

END 
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SUBROUTINE WALK(F, TAIL, NTAIL, HEAD, NEXT, PRED) 

C THIS SUBROUITNE WALKS DOWN THE LINKED LIST AND PLACES THE HEAD (NTAIL) 

C IN THE CORRECT POSITION IN THE HEAD ARRAY. (DECREASING ORDER) 

C INPUT; F, TAIL, NTAIL ARE USED AS ABOVE 

C OUTPUT: ALTERS POINTERS IN NEXT AND PRED ARRAYS TO PLACE HEAD(NTAIL) 

C AFTER HEAD(F) AND BEFORE HEAD (TAIL) IN DOUBLE LINK LIST 

IMPLICIT INTEGER (A-Z) 

DIMENSION HEAD(IOOO), NEXT(IOOO), PRED(IOOO) 

I = F 

NEXT I = NEXT(I) 

HNTAIL = HEAD(NTAIL) 

1 IF(HNTAIL.GE.HEAD(NEXTI )) THEN 
NEXT(I) = NTAIL 
PRED( NTAIL) = I 
NNTAIL = NEXT(NTAIL) 

NEXT (NTAIL) = NEXT I 
PRED (NEXT I ) = NTAIL 
NEXT(TAIL) = NNTAIL 
PRED(NNTAIL) = TAIL 
RETURN 
END IF 

I = NEXT I 
NEXT I = NEXT(I) 

GOTO 1 

END 
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SUBROUTINE CANCEL (START, HEAD, NEXT, PRED, NELEMS) 

C GIVEN THE SORTED RAW STRING, THIS SUBROUTINE REDUCES THE STRING TO 
C CANONICAL FORM BY ADJUSTING THE LINK LIST TO SKIP ANY PAIRS OF 
C ADJACENT LIKE ELEMENTS IN THE STRING. ONCE TWO ELEMENTS ARE 
C ELIMINATED FROM THE STRING, THE SUBSTRINGS WHICH WERE ON THEIR LEFT 
C AND RIGHT CONCATENATE TO FORM A NEW STRING AND HENCE A NEW PAIR OF 
C ADJACENT ELEMENTS WHICH MUST ALSO BE CHECKED FOR EQUALITY. IN THE 
C EVENT A BETA SUBSTRING BECOMES ANNHILATED FROM THE STRING, THE TWO 
C ADJACENT ALPHA STRINGS CONCATENATE AND ARE AGAIN SORTED FROM SMALLEST 
C TO LARGEST IN ABSOLUTE VALUE AS A SINGLE SUBSTRING. 

C INPUT: SORTED RAW STRING 

C OUTPUT: CANONICAL FORM OF THE RAW STRING 

IMPLICIT INTEGER (A-Z) 

DIMENSION HEAD(IOOO), NEXT(IOOO), PRED(IOOO) 

START = 1 
START 1 = START 
START2 = START 
BETA = START 
PTR1 = START 1 
PTR2 = START 1 
50 PTR3 = NEXT(PTR2) 

IF (HEAD(PTR3).LT.O) THEN 
START 1 = PTR3 
START2 = START 1 

ELSE 

BETA = PTR3 
END IF 

10 IF (PTR3.NE.0) THEN 

CALL CHEKER (HEAD, START1, START2, PTR1, PTR2, PTR3, BETA) 

IF (HEAD(PTR2).EQ.HEAD(PTR3)) THEN 
NEXT(PTRI) = NEXT(PTR3) 

PRED(NEXT(PTR3)) = PTR1 
PTR2 = NEXT(PTR3) 

IF (PTR2.EQ.0) GO TO 60 
CHECK1 = HEAD(PTR2) * HEAD(PTR3) 

CHECK2 = HEAD(PTRI) * HEAD(PTR2) 

IF (CHECK1.LT.0.AND.CHECK2.GT.0.AND.HEAD(PTR2).LT.0)THEN 
PTR2 = PTR1 
PTR1 = PRED(PTRI) 

PTR3 = NEXT(PTR2) 

START 1 = START2 

6 I F ( START 1 . NE . 1 . AND . HEAD( PRED ( START 1 ) ) . LT . 0 ) THEN 
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START1 = PRED(STARTI) 

00 TO 6 
END IF 

START2 = PTR3 

CALL MER(X (HEAD, NEXT, START1, START2, PRED) 

PTR1 = PRED ( START 1) 

PTR2 = PTR1 
00 TO 50 

ELSEIF (CHECK1.LT.O.AND.CHECK2.EQ,O.AND.HEAD(PTR2).LT.O) 
+ THEN 

PTR2 = PTR1 
PTR3 = NEXT(PTRI) 

START 1 = PTR3 
START2 = START 1 
00 TO 10 

ELSEIF (CHECK1 .OT.O.AND.CHECK2.0T.O.AND.HEAD(PTR2) .OT.O) 
+ THEN 

IF (BETA.EQ.2) THEN 
BETA = 1 

ELSEIF (BETA.EQ.1) THEN 
00 TO 9 

ELSE 

IF(STARTI.NE.I) BETA = PRED(STARTI) 

11 IF(HEAD(PRED(BETA)).OT.O) THEN 

BETA = PRED(BETA) 

00 TO 11 
END IF 
END IF 

9 PTR1 = PRED (BETA) 

PTR2 = BETA 
PTR3 = NEXT(PTR2) 

00 TO 10 

ELSEIF (CHECK1.LT.0. AND. CHECK2. OT.O. AND. HEAD(PTR2). OT.O) 
+ THEN 

PTR2 = PTR1 
PTR1 = PRED(PTRI) 

IF(HEAD(PTR1).LT.0)THEN 
START 1 = PTR1 

5 IF (HEAD(PRED(START1)).LT.0)THEN 

START 1 = PRED( START 1) 

00 TO 5 

END IF 

END IF 

PTR3 = NEXT(PTR2) 

00 TO 10 

ELSEIF ( CHECK 1 . LT . 0 . AND . CHECK2 . LT . 0 . AND . HEAD ( PTR2 ) . LT . 0 ) 
+ THEN 

IF(STARTI.EQ.I) THEN 
START 1 = PTR2 
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START2 = PTR2 
PTR3 = NEXT(PTR2) 

ELSE 

START 1 = START2 
PTR3 = NEXT(PTR2) 

END IF 
QO TO 10 

ELSE 

START 1 = START2 
PTR3 = NEXT(PTR2) 

END IF 

ELSE 

PTR1 = PTR2 
PTR2 = PTR3 
PTR3 = NEXT(PTR2) 

END IF 
00 TO 10 
END IF 

60 CALL COUNTR (START, NEXT, HEAD, NELEMS) 

RETURN 

END 
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SUBROUTINE CHEKER(HEAD, START 1 , START2, PTR1, PTR2, PTR3, BETA) 

C AS POINTERS MOVE ALONG THE CHARACTER STRING DURING CANCELLATION THIS 
C SUBROUTINE CHECKS TO DETERMINE WHETHER OR NOT THE END OF ONE 
C SUBSTRING IS REACHED AND A NEW ONE BEGINS. DEPENDING ON THE OUTCOME 
C OF THIS CHECK, THE POINTERS USED TO IDENTIFY THE BEGINNING OF THE 
C TWO LATEST ALPHA STRINGS AND THE LATEST BETA STRING ARE UPDATED 

INPUT: HEADO ARRAY AND PRESENT POINTER LOCATIONS FROM 'CANCEL' 

SUBROUTINE 

C OUTPUT: MODIFIED INDICES FOR LOCATION OF ALPHA AND BETA STRINGS 

INTEGER HEAD(IOOO), START 1 , START2, PTR1, PTR2, PTR3, CHECK, BETA 
CHECK = HEAD(PTR3) * HEAD(PTR2) 

IF (CHECK. GT.O) THEN 

IF (HEAD(PTR2).LT.O. AND. HEAD(PTRI). GT.O) START2 = PTR2 
IF (HEAD(PTR2). GT.O. AND. HEAD(PTRl).LT.O) BETA = PTR2 
ELSE IF (CHECK. LE.O) THEN 

IF (HEAD(PTR2). GT.O. AND. HEAD(PTRI).LT.O) THEN 
BETA = PTR2 
START 1 = START2 
START2 = PTR3 

ELSE IF (HEAD(PTR2).GT.O) THEN 
IF (START1 .NE. 1) THEN 
START 1 = START 2 
START2 = PTR3 

ELSE 

START 1 = PTR3 
START2 = START 1 
END IF 

ELSEIF (HEAD(PTR2).LT.O. AND. HEAD(PTRl). GT.O) THEN 
BETA = PTR3 
START 1 = START2 
START2 = PTR2 

ELSEIF (HEAD(PTR2).LT.O.AND.HEAD(PTR1).LT.O) THEN 
BETA = PTR3 
END IF 
END IF 

RETURN 

END 
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SUBROUTINE MERGE(HEAD, NEXT, START 1 , START2, PRED) 

Q ♦**:t!***************************************:tc*********:t:********:M:*****# 

C THIS SUBROUTINE MERGES TWO ORDERED ALPHA STRINGS. THE RESULTING 
C SUBSTRING WILL BE INCREASING IN ABSOLUTE VALUE. 

INPUT: HEADO, NEXT() AND PRED() ARRAYS ALONG WITH INDICES START1 

AND START2 INDICATING BEGINNING OF THE ALPHA STRINGS TO BE 
MERGED 

OUTPUT: MODIFIED STRING CONTAINING THE MERGED ALPHA STRINGS WHICH 

WERE A RESULT OF ANNHILATION OF A BETA STRING 

INTEGER HEAD(IOOO), NEXT(IOOO), START1, START2, Pi, P2, TAIL, 

+ PRED (1000) 

LOGICAL GO 

IF (START1 ,EQ. 1) GO TO 2 

IF(HEAD(PRED(START1)).EQ.O) START1 = NEXT(PRED(START1 ) ) 

I F(NEXT(PRED(START1 ) ) . NE . START 1 )THEN 
START1 = NEXT(PRED(START1)) 

7 IF (START1 .EQ. 1) GO TO 2 
I F ( HEAD ( START 1 ). GT . 0 )THEN 

START1 = NEXT(STARTI) 

GO TO 7 
END IF 

8 IF(HEAD(PRED(START1)).LT.0)THEN 
START 1 = PRED( START 1) 

GO TO 8 

END IF 
END IF 

2 CALL SET(HEAD, NEXT, START 1 , START2, PI, P2, TAIL, PRED) 

1 IF(QO(P1, P2, START2, HEAD, NEXT)) THEN 

CALL MERGER(HEAD, NEXT, PI, P2, TAIL, PRED) 

GOTO 1 
END IF 

CALL FXTAIL(HEAD, NEXT, PI, P2, TAIL, START2, PRED) 

RETURN 
END 
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SUBROUTINE SET(HEAD, NEXT, START1, START2, PI, P2, TAIL, PRED) 

Q *dH**t******************************************************Hli****ilPlii*t** 



C SETS THE STARTING POINTER OF THE NEW STRING AND THE TWO POINTERS 
C OF THE STRINGS TO BE MERGED 

C INPUT: POINTERS INTO TWO ALPHA STRINGS 

C OUTPUT: POINTER INTO BEGINNING OF A NEW SORTED STRING 

INTEGER HEAD(IOOO), NEXT(IOOO), START 1 , START2, PI, P2,TAIL, 

+ PRED (1000) 

PI = START 1 
P2 = START2 

IF(HEAD(P1).GT.HEAD(P2)) THEN 
TAIL = PI 
START 1 = PI 
PI = NEXT(PI) 

ELSE 

TAIL = P2 
START 1 = P2 
NEXT(PRED(P1)) = P2 
PRED(P2) = PRED(PI) 

P2 = NEXT(P2) 

END IF 

RETURN 

END 
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LOGICAL FUNCTION QO(P1, P2, START2, HEAD, NEXT) 

C LOGICAL FUNCTION TO DETERMINE WHEN THE END OF EITHER LIST IS REACHED 

C INPUT: POINTERS INTO TWO ALPHA SUBSTRINGS 

OUTPUT: LOGICAL VARIABLE WHICH IS TRUE IF EITHER POINTER IS AT THE 
END OF A STRING, AND FALSE OTHERWISE 

INTEGER HEAD(IOOO), NEXT(IOOO), PI, P2, START2 

GO = .TRUE. 

IF(P1 .EQ.START2)THEN 
GO = .FALSE. 

RETURN 
END IF 

IF(NEXT(P2).EQ.0)THEN 
GO = .FALSE. 

RETURN 
END IF 

IF(HEAD(P2).GT.0)THEN 
GO = .FALSE. 

RETURN 
END IF 

RETURN 
END 
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SUBROUTINE MERGER (HEAD, NEXT, PI, P2, TAIL, PRED) 

C THIS SUBROUTINE COMPARES TWO ELEMENTS - ONE IN EACH ALPHA SUBSTRING 
C AND PLACES THE ONE THAT IS THE SMALLEST ON THE LIST REPRESENTING THE 
C NEW SORTED LIST 

C INPUT: HEAD, NEXT, PRED, AND POINTERS INTO THE HEAD ARRAY 

OUTPUT: HEAD, NEXT, PRED, AND POINTERS FOR THE NEXT TWO ELEMENTS 
TO BE CHECKED 

INTEGER HEAD(IOOO), NEXT(IOOO), PI, P2, TAIL , PRED(IOOO) 

IF(HEAD(P1).GT.HEAD(P2)) THEN 

CALL TACON1(NEXT, PI, TAIL, PRED) 

ELSE 

CALL TACON2(NEXT, P2, TAIL, PRED) 

END IF 

RETURN 
END 



SUBROUTINE TACON1 (NEXT, PI, TAIL, PRED) 

C ADDS THE ELEMENT TO WHICH PI POINTS TO THE TAIL OF THE NEW STRING 



c 


INPUT: 


NEXT, 


PRED, 


PI, 


TAIL 


c 


OUTPUT: 


NEXT, 


PRED, 


PI, 


TAIL 



INTEGER NEXT(IOOO), PI, TAIL, PRED(IOOO) 

NEXT(TAIL) = PI 
PRED(PI) = TAIL 
TAIL = PI 
PI = NEXT(PI) 

RETURN 

END 
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SUBROUTINE TACON2 (NEXT, P2, TAIL, PRED) 

C ADDS THE ELEMENT TO WHICH P2 POINTS TO THE TAIL OF THE NEW STRING 



c 


INPUT: 


NEXT, 


PRED, 


P2, 


TAIL 


c 


OUTPUT: 


NEXT, 


PRED, 


P2, 


TAIL 



INTECXR NEXT (1000), P2, TAIL, PRED(IOOO) 

NEXT(TAIL) = P2 
PRED(P2) = TAIL 
TAIL = P2 
P2 = NEXT(P2) 

RETURN 

END 
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SUBROUTINE FXTAIL(HEAD, NEXT, PI, P2, TAIL, START2, PRED) 

Q ********♦*♦****♦***♦*:*:♦♦♦*♦♦♦♦*♦*♦*♦**♦♦*♦*♦♦*♦♦♦**♦♦*♦♦♦*♦*♦*♦♦♦♦♦♦*♦ 

C ATTACHES THE END OF THE LONGER SORTED STRING TO THE END OF THE MERCED 
C STRING 

C INPUT: PI, P2, HEAD, NEXT, PRED, TAIL 

C OUTPUT: PI, P2, HEAD, NEXT, PRED, TAIL, START2 

INTEGER HEAD(IOOO), NEXT(IOOO), PI, P2, TAIL, START2, P11, 

+ PRED(IOOO) 

IF(P1 .EQ.START2)THEN 

CALL TACON2(NEXT, P2, TAIL, PRED) 

RETURN 
END IF 

IF(NEXT(P2).EQ.0)THEN 

CALL TAOON1(NEXT, PI, TAIL, PRED) 

P11 = TAIL 

2 IF(P1.NE,START2)THEN 

P11 = PI 
PI = NEXT(PI) 

GOTO 2 
END IF 

NEXT(PII) = P2 
PRED(P2) = P11 
RETURN 
END IF 

IF(HEAD(P2).GT.0)THEN 

CALL TACONKNEXT, PI, TAIL, PRED) 

P11 = TAIL 

3 IF(P1 .NE.START2) THEN 

P11 = PI 
PI = NEXT(PI) 

GOTO 3 
ENDIF 

NEXT(PII) = P2 
PRED(P2) = P11 
RETURN 
ENDIF 

RETURN 

END 
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SUBROUTINE CTREX(FL, NL, PL, FR, NR, PR, BSEED, PSEED) 

C THIS SUBROUTINE CHECKS F(C(R(P))) AGAINST K(F(R(P))) FOR EQUALITY. 

C INPUT: FUNDAMENTAL GROUP REPRESENTATIONS FOR BOTH ALGORITHMS 

OUTPUT: MESSAGE INDICATING THE OCCURRENCE OF UNEQUAL FUNDAMENTAL 

GROUPS 

INTEGER FL(IOOO), NL(IOOO), PL(IOOO), FR(IOOO), NR(IOOO), PR(IOOO) 

IF(NL(1).EQ.O.AND.NR(1).EQ.O) GO TO 11 

I = NL(1) 

J = NR(1) 

10 IF (I .EQ.O.AND.J.EQ.O) GOTO 11 

IF (FL( I ).EQ.O.AND.FR(J).EQ.O) GO TO 11 

IF (FL( I ).EQ.FR(J)) THEN 
I = NL(I) 

J = NR(J) 

GOTO 10 

ELSE 

PR I NT* , ' STR I NGS UNEQUAL ' 

PRINT*,' ' 

IF(NL(1).EQ.0)THEN 

PRINT*, 'THE GROUP IS EMPTY AFTER ALG1 SEQUENCE' 

ELSE 

PRINT*, 'THE FUNDAMENTAL GROUP AFTER ALG1 IS:' 

CALL PRINTS(FL,NL) 

END IF 

IF(NR(1).EQ.O) THEN 

PRINT*, 'THE GROUP IS EMPTY AFTER ALG2 SEQUENCE' 

ELSE 

PRINT*, 'THE FUNDAMENTAL GROUP AFTER ALG2 IS:' 

CALL PRINTS(FR,NR) 

END IF 
PRINT*,' ' 

PRINT*, 'BSEED AND PSEED ARE' , BSEED, PSEED 
END IF 

11 RETURN 
END 
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SUBROUTINE PR I NTS (HEAD, NEXT) 

C THIS SUBROUTINE CAN BE USED TO PRINT OUT ANY STRING THAT IS STORED IN 
C LINK LIST FORM 

C INPUT; LINK LIST ARRAYS 

C OUTPUT: HORIZONTAL STRING OF THE ELEMENTS ACCORDING TO THE NEXT() 

INTEGER HEAD(IOOO), NEXT(IOOO), STRING(IOOO) 

NSTART = NEXT(1) 

I = 0 

98 IF(NSTART.NE.O) THEN 
1 = 1 + 1 

STRING(I) = HEAD(NSTART) 

NSTART = NEXT (NSTART) 

GOTO 98 
END IF 

PRINT 111, (STRING(J), J = 1, I) 

111 FORMATC ’,2014) 

PRINT*,' ’ 

RETURN 

END 
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C ALL SUBROUTINES BELOW THIS POINT WERE USED AS DEBUG TOOLS. THEY 
C GENERATE CRUDE GRAPHS WHICH PLOT THE BOARDS AND THE PATHS. 

C THE VALUES AQUIRED BY THE RAWSTRING , MERGE, AND MANY OTHER 
C SUBROUTINES WERE CHECKED USING THESE ROUTINES. 

Q ♦♦********:^**********************************:^***********:(c**********:t: 

SUBROUTINE GRAPH (BX, BY, X, Y, NOBS, NUMPTS) 

Q :(c4:4::4c:tc4:4::tc:tc4::(c:tc:|c4:4::tc:|c4::tc:tc4:4::4c:tc:4c:4c:tc:tc4::4c:|c:|c:)::|c:tc:(c:t;:4c:tc4E4::4c:tc:|c4:4::tc:|c4::4c:tc:|c4;4;:|c:tc4:4::(c:tc4:4::tc:tc4:4::4c:tc4::(c 

REAL*8 X(1000), Y(1000), BX(IOOO), BY(IOOO) 

CHARACTER*! MATR(53, 105) 

CALL GRAFPA(X, Y, NUMPTS, MATR) 

CALL GRAFOB(BX, BY, NOBS, MATR) 

DO 1 I = 1,53 

1 PRINT 111, (MATR(I,J), J=1,105) 

111 FORMATC ',105A1) 

RETURN 

END 
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SUBROUTINE GRAPPA (X, Y, NUMPTS, MATR) 

REAL*8 X(1000), Y(1000) 

CHARACTER*! MATR (53, 105), CH 

CX) 10 L = 2,52 

DO 10 M = 2,104 
10 MATR(L,M) = ’ ’ 

DO 1 I = 1, 105 

MATR(1,I) = '*’ 

1 MATR(53, I ) = ’*’ 

DO 4 I = 1, 53 

MATR(I,1) = ’*’ 

4 MATR(I ,105) = •*’ 

XK = X(1) 

YK = Y(1) 

CALL COOR(XK, YK, IX, lY) 

CALL CHARPA(1, CH) 

MATR( lY, IX) = CH 

DO 2 K = 2, NUMPTS 

2 CALL FILL(K, X, Y, MATR) 

DO 3 K = 2, NUMPTS 
XK = X(K) 

YK = Y(K) 

CALL COOR(XK, YK, IX, lY) 

CALL CHARPA(K, CH) 

3 MATR( lY, IX) = CH 

RETURN 

END 



SUBROUTINE COOR(BX, BY, IX, lY) 

REAL*8 BX, BY 

IX = MESHP(BX) * 2 
lY = 54 - MESHP(BY) 

RETURN 

END 
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INTEGER FUNCTION MESHP(X) 
REAL*8 X 

MESHP = ((X+1.)*(51./2.)) + 2. 

RETURN 

END 



SUBROUTINE CHARPA(K, CH) 

CHARACTER*! CH 

IF(K.EQ.I) THEN 
CH = 'A' 

RETURN 

END IF 

IF(K.EQ.2) THEN 
CH = 'B' 

RETURN 

END IF 

IF(K.EQ.3) THEN 
CH = 'C 
RETURN 

END IF 

IF(K.EQ.4) THEN 
CH = 'D' 

RETURN 

END IF 

IF(K.EQ.5) THEN 
CH = 'E' 

RETURN 

END IF 

IF(K.EQ.6) THEN 
CH = 'A' 

RETURN 

END IF 

END 



98 



SUBROUTINE FILL(K, X, Y, MATR) 

REAL*8 X(1000), Y(1000), U, V, DX, DY, XK, YK, XKM1 , YKM1 
CHARACTER*! MATR (53, 105) 

XK = X(K) 

YK = Y(K) 

XKM1 = X(K-1) 

YKM1 = Y(K-1) 

CALL COOR(XK, YK, KX, KY) 

CALL OOOR(XKM1, YKM1 , KXM1 , KYM1) 

L = IABS(KX - KXM1) - 1 
M = IABS(KY - KYM1) - 1 
IF (M.GT.L) L = M 
DX = (XK - XKM1)/(L+1) 

DY = (YK - YKM1)/(L+1) 

DO 1 J = 1,L 

U = XKM1 + J*DX 
V = YKM1 + J*DY 
CALL COOR(U,V, lU, IV) 

1 MATR(IV, lU) = ’+’ 

RETURN 

END 



SUBROUTINE GRAFOB (X, Y, NOBS, MATR) 

REAL*8 X(1000), Y(1000) 

CHARACTER*! MATR (53, !05), CH 

MATR(27, 53) = ’*' 

DO 2 K = !,NOBS 
XK = X(K) 

YK = Y(K) 

CALL COOR(XK, YK, IX, lY) 

CALL CHAR(K, CH) 

2 MATR( lY, IX) = CH 

RETURN 

END 
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SUBROUTINE CHAR(K, CH) 

CHARACTER*! CH 

IF(K.EQ.I) THEN 
CH = '1' 

RETURN 
END IF 

IF(K.EQ.2) THEN 
CH = ’2’ 

RETURN 
END IF 

IF(K.EQ.3) THEN 
CH = *3' 

RETURN 
END IF 

IF(K.EQ.4) THEN 
CH = ’4’ 

RETURN 
END IF 

IF(K.EQ.S) THEN 
CH = ’5* 

RETURN 
END IF 

IF(K.EQ.6) THEN 
CH = '6' 

RETURN 
END IF 

IF(K.EQ.7) THEN 
CH = '7' 

RETURN 
END IF 

IF(K.EQ.S) THEN 
CH = '8' 

RETURN 
END IF 

IF(K.EQ.9) THEN 
CH = ’9’ 

RETURN 
END IF 

IF(K.EQ.IO) THEN 
CH = 

RETURN 
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END IF 



IF(K.EQ. 11) THEN 
CH = 

RETURN 
END IF 

IF(K.EQ.12) THEN 
CH = ■$’ 
RETURN 
END IF 

IF(K.EQ.13) THEN 
CH = '%' 
RETURN 
END IF 

IF(K.EQ.14) THEN 
CH = 

RETURN 
END IF 

IF(K.EQ.15) THEN 
CH = ‘4:' 
RETURN 
END IF 

IF(K.EQ.16) THEN 
CH = ■/' 
RETURN 
END IF 

IF(K.EQ. 17) THEN 
CH = '\' 
RETURN 
END IF 

IF(K.EQ.IS) THEN 
CH = ’<’ 
RETURN 
END IF 

IF(K.EQ.19) THEN 
CH = ’>' 
RETURN 
END IF 

IF(K.EQ.20) THEN 
CH = '?■ 
RETURN 
END IF 
END 
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APPENDIX B. A SMALL EXAMPLE 



Given the following five classes (pictured in Figure B.l) 
as candidates for the class containing the true shortest path: 
We formulate a lower bound for each class. 

1 . 

2 • ct2lB^ 

3 • Q|32 

4 • Q[2l^j^^2^2 

5 . a^B2 




Figure B.l The Five Classes Associated With Two Obstacles 
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Consider the class a|I32 pictured below in Figure B.2. 




Figure B.2 The Class o:^B2 



To formulate the lower bound for each class we first 
shrink the obstacles down to a point (Figure B.3). In order 
to cross a 2 we develop a cone of directions (Figure B.3) in 
which the first segment of the bounding polygonal path must 




Figure B.3 The Cone of Directions Associated With qj 
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This cone is then intersected with the one associated with an 



132 crossing pictured in Figure B.4. 




Figure B.4 Cone Associated With B 2 



The resulting cone of directions is denoted by the double arc 
in Figure B.5. 




Figure B.5 The Intersected Cone 
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since we have reached the end of the class name we check 



to see if z lies in this cone of directions. If z lies in the 
cone we go directly to z. Otherwise, we determine which side 
of the cone z lies so that we can find the peg on which the 
path bends. 

In this example, z does not lie in the final cone of 
directions so the path bends at the point bj . Thus the 
bounding path is made of two line segments. One from a to b 2 
and the other from b 2 to z (Figure B.6). 



Figure B.6 The Bounding Path for the Class 

A similar procedure is done on each of the remaining four 
classes yielding the following lower bounds: 




2 



1 . 0202 



1.5231 




1 . 5303 




2 . 2885 



4 . O2B2 ^2^2 



3.2885 



5 . o 2 B 2 



2.7008 
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This list is then put into increasing order: 




1.5231 




1.5303 




2.2885 



4 . 



2.7008 



5 . ct2B^ct2l^2 



3 . 2885 



Now the first class on this ordered list is chosen. To 
find the true shortest path we expand the obstacles back to 
their original shape. If the path bends on a peg, then a 
tangent line is constructed from a start point to the correct 
side of the obstacle (Figure B.7). 



In the class a|B 2 / must draw a tangent line from a to 
b 2 initially. The point of tangency on B 2 can be found by 
using the formulas developed in Chapter V of this paper. We 
then travel around the obstacle until a line tangent to B 2 can 
be drawn to z (Figure B.8). 




Figure B.7 A Tangent Line From a to B 2 
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Throughout this process the lengths of the line segments 
and the arc lengths representing the distance traveled around 
the obstacles is summed. The final sum is the length of the 
true shortest path for the respective class. 

If the shortest path in a class is smaller than the lower 
bound on the next class to be considered, then the search 
stops. In this example the search stops after the first class 
is checked since the length of o^ct 2 is 1.8793 and the length 
of a 2 B, is 1.9032. 
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